Merge pull request #69 from vavrines/dev

Dev

Project.toml

@@ -1,7 +1,7 @@
 name = "KitBase"
 uuid = "86eed249-3a28-466f-8d3a-596821e1af9a"
 authors = ["Tianbai Xiao <[email protected]>"]
-version = "0.9.0"
+version = "0.9.1"
 
 [deps]
 BSON = "fbb218c0-5317-5bc6-957e-2ee96dd4b1f0"

README.md

@@ -5,7 +5,7 @@
 [![](https://img.shields.io/badge/docs-stable-blue.svg)](https://xiaotianbai.com/Kinetic.jl/stable/)
 [![codecov](https://img.shields.io/codecov/c/github/vavrines/KitBase.jl)](https://codecov.io/gh/vavrines/KitBase.jl)
 [![deps](https://juliahub.com/docs/KitBase/deps.svg)](https://juliahub.com/ui/Packages/KitBase/YOFTS?t=2)
-[![GitHub commits since tagged version](https://img.shields.io/github/commits-since/vavrines/KitBase.jl/v0.8.0.svg?style=social&logo=github)](https://github.com/vavrines/KitBase.jl)
+[![GitHub commits since tagged version](https://img.shields.io/github/commits-since/vavrines/KitBase.jl/v0.9.0.svg?style=social&logo=github)](https://github.com/vavrines/KitBase.jl)
 
 This lightweight module provides basic physical formulations and numerical methods for [Kinetic.jl](https://github.com/vavrines/Kinetic.jl) ecosystem.
 The finite volume method (FVM) is employed to perform 1-3 dimensional numerical simulations on CPUs and GPUs.

src/Flux/flux_boundary.jl

@@ -15,14 +15,14 @@ Maxwell's diffusive boundary flux
 * ``rot``: rotation indicator (1/-1)
 """
 function flux_boundary_maxwell!(
-    fw::AV{T1},
-    bc::AV{T2},
-    w::AV{T3},
+    fw::AV,
+    bc::AV,
+    w::AV,
     inK,
     γ,
     dt,
     rot, # 1 / -1
-) where {T1<:FN,T2<:Real,T3<:Real}
+)
 
     @assert length(bc) == 3
 
@@ -48,18 +48,18 @@ $(SIGNATURES)
 1D2F1V
 """
 function flux_boundary_maxwell!(
-    fw::AV{T1},
-    fh::T2,
-    fb::T2,
-    bc::AV{T3},
-    h::T4,
-    b::T4,
-    u::T5,
-    ω::T5,
+    fw::AV,
+    fh::T1,
+    fb::T1,
+    bc::AV,
+    h::T2,
+    b::T2,
+    u::T3,
+    ω::T3,
     inK,
     dt,
     rot = 1,
-) where {T1<:FN,T2<:AV{<:FN},T3<:Real,T4<:AV{<:FN},T5<:AV{<:FN}} # 1D2F1V
+) where {T1<:AV,T2<:AV,T3<:AV}
 
     @assert length(bc) == 3
 
@@ -92,18 +92,67 @@ end
 """
 $(SIGNATURES)
 
-2D Continuum
+Mixture
 """
 function flux_boundary_maxwell!(
-    fw::AV{T1},
-    bc::AV{T2},
-    w::AV{T3},
+    fw::AM,
+    fh::T1,
+    fb::T1,
+    bc::AM,
+    h::T2,
+    b::T2,
+    u::T3,
+    ω::T3,
     inK,
-    γ,
     dt,
-    len,
-    rot,
-) where {T1<:Real,T2<:Real,T3<:Real}
+    rot = 1,
+) where {T1<:AM,T2<:AM,T3<:AM}
+
+    @assert size(bc, 1) == 3
+
+    δ = heaviside.(u .* rot)
+    SF = [sum(ω[:, j] .* u[:, j] .* h[:, j] .* (1.0 .- δ[:, j])) for j in axes(h, 2)]
+    SG = [
+        (bc[end, j] / π)^0.5 * sum(
+            ω[:, j] .* u[:, j] .* exp.(-bc[end, j] .* (u[:, j] .- bc[2, j]) .^ 2) .*
+            δ[:, j],
+        ) for j in axes(h, 2)
+    ]
+
+    prim = zero(bc)
+    for j in axes(prim, 2)
+        prim[:, j] .= [-SF[j] / SG[j]; bc[2:end, j]]
+    end
+
+    H = mixture_maxwellian(u, prim)
+    B = mixture_energy_maxwellian(H, prim, inK)
+
+    hWall = H .* δ .+ h .* (1.0 .- δ)
+    bWall = B .* δ .+ b .* (1.0 .- δ)
+
+    for j in axes(fw, 2)
+        fw[1, j] = discrete_moments(hWall[:, j], u[:, j], ω[:, j], 1) * dt
+        fw[2, j] = discrete_moments(hWall[:, j], u[:, j], ω[:, j], 2) * dt
+        fw[3, j] =
+            (
+                0.5 * discrete_moments(hWall[:, j] .* u[:, j] .^ 2, u[:, j], ω[:, j], 1) +
+                0.5 * discrete_moments(bWall[:, j], u[:, j], ω[:, j], 1)
+            ) * dt
+    end
+
+    @. fh = u * hWall * dt
+    @. fb = u * bWall * dt
+
+    return nothing
+
+end
+
+"""
+$(SIGNATURES)
+
+2D Continuum
+"""
+function flux_boundary_maxwell!(fw::AV, bc::AV, w::AV, inK, γ, dt, len, rot)
 
     @assert length(bc) == 4
 
@@ -129,17 +178,17 @@ $(SIGNATURES)
 2D1F2V
 """
 function flux_boundary_maxwell!(
-    fw::AV{T1},
-    fh::AM{T2},
-    bc::AV{T3},
-    h::AM{T4},
-    u::T5,
-    v::T5,
-    ω::T5,
+    fw::AV,
+    fh::AM,
+    bc::AV,
+    h::AM,
+    u::T,
+    v::T,
+    ω::T,
     dt,
     len,
     rot = 1,
-) where {T1<:FN,T2<:FN,T3<:Real,T4<:FN,T5<:AM{<:FN}}
+) where {T<:AM}
 
     @assert length(bc) == 4
 
@@ -172,10 +221,10 @@ $(SIGNATURES)
 2D2F2V
 """
 function flux_boundary_maxwell!(
-    fw::AV{T1},
+    fw::AV,
     fh::T2,
     fb::T2,
-    bc::AV{T3},
+    bc::AV,
     h::T4,
     b::T4,
     u::T5,
@@ -185,7 +234,7 @@ function flux_boundary_maxwell!(
     dt,
     len,
     rot = 1,
-) where {T1<:FN,T2<:AM{<:FN},T3<:Real,T4<:AM{<:FN},T5<:AM{<:FN}}
+) where {T2<:AM,T4<:AM,T5<:AM}
 
     @assert length(bc) == 4
 
@@ -227,18 +276,18 @@ $(SIGNATURES)
 1F3V
 """
 function flux_boundary_maxwell!(
-    fw::AV{T1},
-    ff::AA{T2,3},
-    bc::AV{T3},
-    f::AA{T4,3},
-    u::T5,
-    v::T5,
-    w::T5,
-    ω::T5,
+    fw::AV,
+    ff::AA{T1,3},
+    bc::AV,
+    f::AA{T2,3},
+    u::T3,
+    v::T3,
+    w::T3,
+    ω::T3,
     dt,
     area,
     rot = 1,
-) where {T1<:FN,T2<:FN,T3<:Real,T4<:FN,T5<:AA{<:FN,3}}
+) where {T1,T2,T3<:AA{T4,3} where {T4}}
 
     @assert length(bc) == 5
 
@@ -277,14 +326,7 @@ Specular reflection boundary flux
 
 1D1F1V
 """
-function flux_boundary_specular!(
-    fw::AV{T1},
-    ff::AV{T2},
-    f::AV{T3},
-    u::T4,
-    ω::T4,
-    dt,
-) where {T1<:Real,T2<:FN,T3<:FN,T4<:AV{<:FN}}
+function flux_boundary_specular!(fw::AV, ff::AV, f::AV, u::T, ω::T, dt) where {T<:AV}
 
     fWall = similar(f)
     for i in eachindex(f)
@@ -307,15 +349,15 @@ $(SIGNATURES)
 1D2F1V
 """
 function flux_boundary_specular!(
-    fw::AV{T1},
+    fw::AV,
     fh::T2,
     fb::T2,
     h::T3,
     b::T3,
     u::T4,
     ω::T4,
     dt,
-) where {T1<:Real,T2<:AV{<:FN},T3<:AV{<:FN},T4<:AV{<:FN}}
+) where {T2<:AV,T3<:AV,T4<:AV}
 
     hWall = similar(h)
     bWall = similar(b)

src/Flux/flux_fluid.jl

@@ -38,7 +38,7 @@ $(SIGNATURES)
 
 HLL flux
 """
-flux_hll!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) = 
+flux_hll!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) =
     flux_hll!(face, ctrL, ctrR, KS.gas, args...)
 
 function flux_hll!(
@@ -113,7 +113,7 @@ Roe's flux with entropy fix
 _P. L. Roe, Approximate Riemann Solvers, Parameter Vectors and Difference Schemes, Journal of Computational Physics, 43, pp. 357-372._
 (_cf. http://cfdbooks.com/cfdcodes.html_)
 """
-flux_roe!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) = 
+flux_roe!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) =
     flux_roe!(face, ctrL, ctrR, KS.gas, args...)
 
 """

src/Flux/flux_gks.jl

@@ -7,7 +7,7 @@ $(SIGNATURES)
 
 Gas kinetic flux
 """
-flux_gks!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) = 
+flux_gks!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) =
     flux_gks!(face, ctrL, ctrR, KS.gas, args...)
 
 """
@@ -17,14 +17,7 @@ Scalar
 
 `flux_gks` is called since there is no in-place operation for scalar
 """
-function flux_gks!(
-    face::Interface,
-    ctrL,
-    ctrR,
-    gas::Scalar,
-    p,
-    dt = 1.0,
-)
+function flux_gks!(face::Interface, ctrL, ctrR, gas::Scalar, p, dt = 1.0)
 
     dxL, dxR = p[1:2]
 
@@ -49,14 +42,7 @@ $(SIGNATURES)
 
 Gas
 """
-function flux_gks!(
-    face::Interface,
-    ctrL,
-    ctrR,
-    gas::Gas,
-    p,
-    dt = 1.0,
-)
+function flux_gks!(face::Interface, ctrL, ctrR, gas::Gas, p, dt = 1.0)
 
     dxL, dxR = p[1:2]
 
@@ -87,14 +73,7 @@ $(SIGNATURES)
 
 Mixture
 """
-function flux_gks!(
-    face::Interface,
-    ctrL,
-    ctrR,
-    gas::Mixture,
-    p,
-    dt = 1.0,
-)
+function flux_gks!(face::Interface, ctrL, ctrR, gas::Mixture, p, dt = 1.0)
 
     dxL, dxR = p[1:2]
 
@@ -131,11 +110,6 @@ end
 $(SIGNATURES)
 
 Gas kinetic scalar flux
-
-* @args: conservative scalars and their slopes
-* @args: viscosity
-* @args: time step and cell size
-* @return: scalar flux
 """
 function flux_gks(u, μ, dt, a = 0, su = 0.0)
     prim = ifelse(a == 0, conserve_prim(u), conserve_prim(u, a))
@@ -253,10 +227,6 @@ $(SIGNATURES)
 
 Gas kinetic Navier-Stokes flux
 
-* @args: conservative variables and their left/right slopes
-* @args: molecular and thermodynamic parameters
-* @args: time step and cell size
-
 Continuous case
 """
 function flux_gks!(

src/Flux/flux_kcu.jl

@@ -3,7 +3,7 @@ $(SIGNATURES)
 
 Kinetic central-upwind (KCU) flux
 """
-flux_kcu!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) = 
+flux_kcu!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) =
     flux_kcu!(face, ctrL, ctrR, KS.gas, KS.vs, args...)
 
 """

src/Flux/flux_kfvs.jl

@@ -3,7 +3,7 @@ $(SIGNATURES)
 
 Kinetic flux vector splitting (KFVS) flux
 """
-flux_kfvs!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) = 
+flux_kfvs!(KS::AbstractSolverSet, face, ctrL, ctrR, args...) =
     flux_kfvs!(face, ctrL, ctrR, KS.gas, KS.vs, args...)
 
 function flux_kfvs!(
@@ -225,10 +225,6 @@ $(SIGNATURES)
 
 Kinetic flux vector splitting (KFVS) flux
 
-- @args: particle distribution functions and their left/right slopes
-- @args: particle velocity quadrature points and weights
-- @args: time step and cell size
-
 1F1V for pure DOM
 """
 function flux_kfvs!(

src/Math/math_general.jl

@@ -56,9 +56,9 @@ $(SIGNATURES)
 
 Gauss Legendre integral for fast spectral method
 
-* @args: number of quadrature points N, integral range [a, b]
-* @args: quadrature points x & weights w
-
+# Arguments
+* ``N``: number of quadrature points, 
+* ``a, b``: integral range
 """
 function lgwt(N::Integer, a::Real, b::Real)
     x = zeros(N)