kernels for boundary conditions are revised to use in boundary_condition.py. It points out that currently kernels in y direction for solving PDE and BC are diffrent. It needs double check in future.

Author: saeb-faraji-gargari

source/package/soli3d/boundary_condition.py

@@ -8,15 +8,20 @@
 
 # import numpy as np
 # from numpy.typing import NDArray
-from .kernels import kernel_i_jm1, kernel_i_jp1, kernel_im1_j, kernel_ip1_j
+from .kernels import (
+    kernel_i_jm1_boundary,
+    kernel_i_jp1_boundary,
+    kernel_im1_j_boundary,
+    kernel_ip1_j_boundary,
+)
 
 
 def boundary_set(
     phi: Any,
     boundary_loc: Any,
     boundary_type: Any,
-    Dirichlet_boundary_value: Any,
-    Neumann_boundary_value: Any,
+    dirichlet_boundary_value: Any,
+    neumann_boundary_value: Any,
     dx: float,
     dz: float,
 ) -> Any:
@@ -45,11 +50,11 @@ def boundary_set(
     # boundary types are defined in "io_data_process.py/default_boundary_type function"
 
     # ---------------------------------------------------------------------------------------------
-    # phi = lfr.where(
-    #     ((boundary_loc == 1) & (boundary_type == 0)),
-    #     Dirichlet_boundary_value,
-    #     phi,
-    # )
+    phi = lfr.where(
+        ((boundary_loc == 1) & (boundary_type == 0)),
+        dirichlet_boundary_value,
+        phi,
+    )
 
     # # kernel_im1_j   i-1, j
     # kernel_im1_j: NDArray[np.uint8] = np.array(
@@ -121,13 +126,13 @@ def boundary_set(
             (boundary_loc == 1)
             & ((boundary_type == 1) | (boundary_type == 5) | (boundary_type == 6))
         ),
-        lfr.focal_sum(phi, kernel_ip1_j) - (dx * Neumann_boundary_value),
+        lfr.focal_sum(phi, kernel_ip1_j_boundary) - (dx * neumann_boundary_value),
         phi,
     )
 
     phi = lfr.where(
         ((boundary_loc == 1) & (boundary_type == 2)),
-        lfr.focal_sum(phi, kernel_i_jp1) - (dz * Neumann_boundary_value),
+        lfr.focal_sum(phi, kernel_i_jp1_boundary) - (dz * neumann_boundary_value),
         phi,
     )
 
@@ -136,13 +141,13 @@ def boundary_set(
             (boundary_loc == 1)
             & ((boundary_type == 3) | (boundary_type == 7) | (boundary_type == 8))
         ),
-        lfr.focal_sum(phi, kernel_im1_j) + (dx * Neumann_boundary_value),
+        lfr.focal_sum(phi, kernel_im1_j_boundary) + (dx * neumann_boundary_value),
         phi,
     )
 
     phi = lfr.where(
         ((boundary_loc == 1) & (boundary_type == 4)),
-        lfr.focal_sum(phi, kernel_i_jm1) + (dz * Neumann_boundary_value),
+        lfr.focal_sum(phi, kernel_i_jm1_boundary) + (dz * neumann_boundary_value),
         phi,
     )
 

source/package/soli3d/kernels.py

@@ -108,3 +108,47 @@
 #     ],
 #     dtype=np.uint8,
 # )
+
+# ---- NOTE kernel_i_jm1_boundary and kernel_i_jp1_boundary
+# are different kernel_i_jm1 and kernel_i_jp1 -------------
+
+
+# kernel_im1_j   i-1, j
+kernel_im1_j_boundary: NDArray[np.uint8] = np.array(
+    [
+        [0, 0, 0],
+        [1, 0, 0],
+        [0, 0, 0],
+    ],
+    dtype=np.uint8,
+)
+
+# kernel_i_jm1   i, j-1    # Check it. It is changed compared to advection-diffusion  test model.
+kernel_i_jm1_boundary: NDArray[np.uint8] = np.array(
+    [
+        [0, 0, 0],
+        [0, 0, 0],
+        [0, 1, 0],
+    ],
+    dtype=np.uint8,
+)
+
+# kernel_ip1_j   i+1, j
+kernel_ip1_j_boundary: NDArray[np.uint8] = np.array(
+    [
+        [0, 0, 0],
+        [0, 0, 1],
+        [0, 0, 0],
+    ],
+    dtype=np.uint8,
+)
+
+# kernel_i_jp1   i, j+1     # Check it. It is changed compared to advection-diffusion  test model.
+kernel_i_jp1_boundary: NDArray[np.uint8] = np.array(
+    [
+        [0, 1, 0],
+        [0, 0, 0],
+        [0, 0, 0],
+    ],
+    dtype=np.uint8,
+)