1535 particle bed modificaiton and alternative packing algorithms

docs/documentation/case.md

@@ -313,7 +313,7 @@ This is enabled by adding ``'elliptic_smoothing': "T",`` and ``'elliptic_smoothi
 | ---:                 | :----:  | :---                |
 | `num_ibs`            | Integer | Number of immersed boundary patches |
 | `num_stl_models`     | Integer | Number of STL/OBJ model entries in the `stl_models` array |
-| `num_particle_beds`  | Integer | Number of particle bed specifications to generate immersed boundary patches from |
+| `num_particle_clouds` | Integer | Number of particle bed specifications to generate immersed boundary patches from |
 | `ib_neighborhood_radius` | Integer | Parameter that controls the neighborhood size for IB detection. |
 | `geometry`           | Integer | Geometry configuration of the patch.|
 | `x[y,z]_centroid`    | Real    | Centroid of the applied geometry in the [x,y,z]-direction. |

docs/module_categories.json

@@ -38,7 +38,7 @@
             "m_compute_cbc",
             "m_boundary_common",
             "m_ibm",
-            "m_particle_bed",
+            "m_particle_cloud",
             "m_igr",
             "m_ib_patches",
             "m_compute_levelset"

examples/2D_mibm_particle_bed/case.py → examples/2D_mibm_particle_cloud/case.py

@@ -85,19 +85,20 @@
             "collision_time": collision_time,
             "ib_coefficient_of_friction": 0.1,
             # Particle bed: 20 free-floating circles placed randomly in region
-            "num_particle_beds": 1,
-            "particle_bed(1)%x_centroid": bed_x,
-            "particle_bed(1)%y_centroid": bed_y,
-            "particle_bed(1)%z_centroid": 0.0,
-            "particle_bed(1)%length_x": bed_lx,
-            "particle_bed(1)%length_y": bed_ly,
-            "particle_bed(1)%length_z": 0.0,
-            "particle_bed(1)%num_particles": 20,
-            "particle_bed(1)%radius": particle_radius,
-            "particle_bed(1)%mass": particle_mass,
-            "particle_bed(1)%min_spacing": particle_min_spacing,
-            "particle_bed(1)%moving_ibm": 2,
-            "particle_bed(1)%seed": 42,
+            "num_particle_clouds": 1,
+            "particle_cloud(1)%x_centroid": bed_x,
+            "particle_cloud(1)%y_centroid": bed_y,
+            "particle_cloud(1)%z_centroid": 0.0,
+            "particle_cloud(1)%length_x": bed_lx,
+            "particle_cloud(1)%length_y": bed_ly,
+            "particle_cloud(1)%length_z": 0.0,
+            "particle_cloud(1)%num_particles": 20,
+            "particle_cloud(1)%radius": particle_radius,
+            "particle_cloud(1)%mass": particle_mass,
+            "particle_cloud(1)%min_spacing": particle_min_spacing,
+            "particle_cloud(1)%moving_ibm": 2,
+            "particle_cloud(1)%seed": 42,
+            "particle_cloud(1)%packing_method": 1,
             # Output
             "format": 1,
             "precision": 2,

src/common/m_constants.fpp

@@ -31,7 +31,7 @@ module m_constants
     !> Fixed capacity of patch_ib (namelist patches + local particle bed subset after reduction)
     integer, parameter  :: num_ib_patches_max_namelist = 54000
     integer, parameter  :: num_local_ibs_max = 2000            !< Maximum number of immersed boundary patches (patch_ib)
-    integer, parameter  :: num_particle_beds_max = 10          !< Maximum number of particle bed patch specifications
+    integer, parameter  :: num_particle_clouds_max = 10        !< Maximum number of particle bed patch specifications
     integer, parameter  :: num_bc_patches_max = 10             !< Maximum number of boundary condition patches
     integer, parameter  :: max_2d_fourier_modes = 10           !< Max Fourier mode index for 2D modal patch (geometry 13)
     integer, parameter  :: max_sph_harm_degree = 5             !< Max degree L for 3D spherical harmonic patch (geometry 14)

src/common/m_derived_types.fpp

@@ -320,7 +320,7 @@ module m_derived_types
         real(wp), dimension(1:3) :: step_angular_vel  !< velocity array used to store intermediate steps in the time_stepper module
     end type ib_patch_parameters
 
-    type particle_bed_parameters
+    type particle_cloud_parameters
         real(wp) :: x_centroid, y_centroid, z_centroid  !< Center of the particle bed region
         real(wp) :: length_x, length_y, length_z  !< Dimensions of the particle bed region
         integer  :: num_particles  !< Number of particles to generate
@@ -329,7 +329,8 @@ module m_derived_types
         real(wp) :: min_spacing  !< Minimum surface-to-surface gap (particle centers are 2*radius + min_spacing apart)
         integer  :: moving_ibm  !< Motion flag: 0=static, 1=moving (forces), 2=forced path
         integer  :: seed  !< Random seed for reproducible placement
-    end type particle_bed_parameters
+        integer  :: packing_method  !< Packing algorithm: 1=rejection sampling
+    end type particle_cloud_parameters
 
     !> Derived type annexing the physical parameters (PP) of the fluids. These include the specific heat ratio function and liquid
     !! stiffness function.

src/simulation/m_checker.fpp

@@ -38,6 +38,10 @@ contains
 
         @:PROHIBIT(ib_state_wrt .and. .not. ib, "ib_state_wrt requires ib to be enabled")
 
+        if (num_particle_clouds > 0) then
+            call s_check_inputs_particle_clouds
+        end if
+
     end subroutine s_check_inputs
 
     !> Checks constraints on compiler options
@@ -106,4 +110,21 @@ contains
 
     end subroutine s_check_inputs_nvidia_uvm
 
+    !> Checks that each active particle cloud has a valid packing_method specified
+    impure subroutine s_check_inputs_particle_clouds
+
+        integer          :: i
+        character(len=5) :: idxStr
+
+        do i = 1, num_particle_clouds
+            call s_int_to_str(i, idxStr)
+            @:PROHIBIT(particle_cloud(i)%packing_method == dflt_int, &
+                       & "particle_cloud("//trim(idxStr)//")%packing_method must be specified (1 = rejection sampling)")
+            @:PROHIBIT(particle_cloud(i)%packing_method /= 1, &
+                       & "particle_cloud("//trim(idxStr) &
+                       & //")%packing_method must be 1 (rejection sampling is the only supported method)")
+        end do
+
+    end subroutine s_check_inputs_particle_clouds
+
 end module m_checker

src/simulation/m_global_parameters.fpp

@@ -263,7 +263,7 @@ module m_global_parameters
     !> @{
     type(ib_patch_parameters), dimension(num_ib_patches_max_namelist) :: patch_ib  !< Immersed boundary patch parameters
     integer, dimension(num_local_ibs_max) :: local_ib_patch_ids  !< lookup table of IBs in the local compute domain
-    type(particle_bed_parameters), dimension(num_particle_beds_max) :: particle_bed  !< Particle bed specifications
+    type(particle_cloud_parameters), dimension(num_particle_clouds_max) :: particle_cloud  !< Particle bed specifications
     integer, allocatable, dimension(:,:,:) :: ib_neighbor_ranks  !< MPI ranks of neighborhood domains, indexed (-N:N,-N:N,-N:N)
     type(ib_airfoil_parameters), dimension(num_ib_airfoils_max) :: ib_airfoil  !< Per-airfoil NACA user inputs (namelist)
     type(ib_airfoil_grid), dimension(num_ib_airfoils_max) :: ib_airfoil_grids  !< Per-airfoil computed surface grids
@@ -689,20 +689,21 @@ contains
             ib_airfoil_grids(i)%Np = 0
         end do
 
-        num_particle_beds = 0
-        do i = 1, num_particle_beds_max
-            particle_bed(i)%x_centroid = 0._wp
-            particle_bed(i)%y_centroid = 0._wp
-            particle_bed(i)%z_centroid = 0._wp
-            particle_bed(i)%length_x = dflt_real
-            particle_bed(i)%length_y = dflt_real
-            particle_bed(i)%length_z = dflt_real
-            particle_bed(i)%num_particles = 0
-            particle_bed(i)%radius = dflt_real
-            particle_bed(i)%mass = dflt_real
-            particle_bed(i)%min_spacing = 0._wp
-            particle_bed(i)%moving_ibm = 0
-            particle_bed(i)%seed = 0
+        num_particle_clouds = 0
+        do i = 1, num_particle_clouds_max
+            particle_cloud(i)%x_centroid = 0._wp
+            particle_cloud(i)%y_centroid = 0._wp
+            particle_cloud(i)%z_centroid = 0._wp
+            particle_cloud(i)%length_x = dflt_real
+            particle_cloud(i)%length_y = dflt_real
+            particle_cloud(i)%length_z = dflt_real
+            particle_cloud(i)%num_particles = 0
+            particle_cloud(i)%radius = dflt_real
+            particle_cloud(i)%mass = dflt_real
+            particle_cloud(i)%min_spacing = 0._wp
+            particle_cloud(i)%moving_ibm = 0
+            particle_cloud(i)%seed = 0
+            particle_cloud(i)%packing_method = dflt_int
         end do
 
         do i = 1, num_ib_patches_max_namelist

src/simulation/m_mpi_proxy.fpp

@@ -74,7 +74,7 @@ contains
             & 'wave_speeds', 'avg_state', 'precision', 'bc_x%beg', 'bc_x%end', &
             & 'bc_y%beg', 'bc_y%end', 'bc_z%beg', 'bc_z%end',  'fd_order',     &
             & 'num_probes', 'num_integrals', 'bubble_model', 'thermal',        &
-            & 'num_source', 'relax_model', 'num_ibs', 'num_particle_beds', 'n_start',    &
+            & 'num_source', 'relax_model', 'num_ibs', 'num_particle_clouds', 'n_start',    &
             & 'num_bc_patches', 'num_igr_iters', 'num_igr_warm_start_iters', &
             & 'adap_dt_max_iters', 'collision_model', 'ib_neighborhood_radius', &
             & 'int_comp' ]
@@ -222,14 +222,15 @@ contains
             #:endfor
         end do
 
-        do i = 1, num_particle_beds
+        do i = 1, num_particle_clouds
             #:for VAR in ['x_centroid', 'y_centroid', 'z_centroid', 'length_x', 'length_y', 'length_z', &
                 & 'radius', 'mass', 'min_spacing']
-                call MPI_BCAST(particle_bed(i)%${VAR}$, 1, mpi_p, 0, MPI_COMM_WORLD, ierr)
+                call MPI_BCAST(particle_cloud(i)%${VAR}$, 1, mpi_p, 0, MPI_COMM_WORLD, ierr)
             #:endfor
-            call MPI_BCAST(particle_bed(i)%num_particles, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
-            call MPI_BCAST(particle_bed(i)%moving_ibm, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
-            call MPI_BCAST(particle_bed(i)%seed, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
+            call MPI_BCAST(particle_cloud(i)%num_particles, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
+            call MPI_BCAST(particle_cloud(i)%moving_ibm, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
+            call MPI_BCAST(particle_cloud(i)%seed, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
+            call MPI_BCAST(particle_cloud(i)%packing_method, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr)
         end do
 
         do j = 1, num_probes_max