From 3426504b845ee798a74f91fe3ff1726fbcdd407a Mon Sep 17 00:00:00 2001 From: Sam Reeve <6740307+streeve@users.noreply.github.com> Date: Wed, 29 May 2024 15:12:54 -0400 Subject: [PATCH 1/2] Ensure profile output works if points are not at exactly zero --- src/CabanaPD_DisplacementProfile.hpp | 36 ++++++++++++++++++++++------ 1 file changed, 29 insertions(+), 7 deletions(-) diff --git a/src/CabanaPD_DisplacementProfile.hpp b/src/CabanaPD_DisplacementProfile.hpp index e0a652ed..39ae06e9 100644 --- a/src/CabanaPD_DisplacementProfile.hpp +++ b/src/CabanaPD_DisplacementProfile.hpp @@ -42,30 +42,52 @@ void createOutputProfile( MPI_Comm comm, const int num_cell, Kokkos::View count( "c", 1 ); auto dims = getDim( profile_dim ); - double dx1 = particles.dx[dims[0]]; - double dx2 = particles.dx[dims[1]]; + double dx1 = particles.dx[dims[0]] / 2.0; + double dx2 = particles.dx[dims[1]] / 2.0; + double center1 = particles.local_mesh_lo[dims[0]] + + particles.global_mesh_ext[dims[0]] / 2.0; + double center2 = particles.local_mesh_lo[dims[1]] + + particles.global_mesh_ext[dims[1]] / 2.0; auto x = particles.sliceReferencePosition(); auto u = particles.sliceDisplacement(); + // Find points closest to the center line. + double center_min1; + double center_min2; + auto find_profile = + KOKKOS_LAMBDA( const int pid, double& min1, double& min2 ) + { + if ( Kokkos::abs( x( pid, dims[0] ) - center1 ) < min1 ) + min1 = x( pid, dims[0] ); + if ( Kokkos::abs( x( pid, dims[1] ) - center2 ) < min2 ) + min2 = x( pid, dims[1] ); + }; + Kokkos::RangePolicy policy( + 0, x.size() ); + Kokkos::parallel_reduce( "displacement_profile", policy, find_profile, + Kokkos::Min( center_min1 ), + Kokkos::Min( center_min2 ) ); + + // Extract points along that line. auto measure_profile = KOKKOS_LAMBDA( const int pid ) { - if ( x( pid, dims[0] ) < dx1 / 2.0 && x( pid, dims[0] ) > -dx1 / 2.0 && - x( pid, dims[1] ) < dx2 / 2.0 && x( pid, dims[1] ) > -dx2 / 2.0 ) + if ( x( pid, dims[0] ) - center_min1 < dx1 && + x( pid, dims[0] ) - center_min1 > -dx1 && + x( pid, dims[1] ) - center_min2 < dx2 && + x( pid, dims[1] ) - center_min2 > -dx2 ) { auto c = Kokkos::atomic_fetch_add( &count( 0 ), 1 ); profile( c, 0 ) = x( pid, profile_dim ); profile( c, 1 ) = user( pid ); } }; - Kokkos::RangePolicy policy( - 0, x.size() ); Kokkos::parallel_for( "displacement_profile", policy, measure_profile ); auto count_host = Kokkos::create_mirror_view_and_copy( Kokkos::HostSpace{}, count ); auto profile_host = Kokkos::create_mirror_view_and_copy( Kokkos::HostSpace{}, profile ); - std::fstream fout; + std::fstream fout; fout.open( file_name, std::ios::app ); for ( int p = 0; p < count_host( 0 ); p++ ) { From 9aa62dd1bc2401f36336ea3e6cd4469895b1e818 Mon Sep 17 00:00:00 2001 From: Sam Reeve <6740307+streeve@users.noreply.github.com> Date: Wed, 29 May 2024 16:18:29 -0400 Subject: [PATCH 2/2] fixup: use global bounds for profile output --- src/CabanaPD_DisplacementProfile.hpp | 4 ++-- src/CabanaPD_Particles.hpp | 6 ++++++ 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/src/CabanaPD_DisplacementProfile.hpp b/src/CabanaPD_DisplacementProfile.hpp index 39ae06e9..bd65efcd 100644 --- a/src/CabanaPD_DisplacementProfile.hpp +++ b/src/CabanaPD_DisplacementProfile.hpp @@ -44,9 +44,9 @@ void createOutputProfile( MPI_Comm comm, const int num_cell, auto dims = getDim( profile_dim ); double dx1 = particles.dx[dims[0]] / 2.0; double dx2 = particles.dx[dims[1]] / 2.0; - double center1 = particles.local_mesh_lo[dims[0]] + + double center1 = particles.global_mesh_lo[dims[0]] + particles.global_mesh_ext[dims[0]] / 2.0; - double center2 = particles.local_mesh_lo[dims[1]] + + double center2 = particles.global_mesh_lo[dims[1]] + particles.global_mesh_ext[dims[1]] / 2.0; auto x = particles.sliceReferencePosition(); diff --git a/src/CabanaPD_Particles.hpp b/src/CabanaPD_Particles.hpp index 02d691c2..d0ea4397 100644 --- a/src/CabanaPD_Particles.hpp +++ b/src/CabanaPD_Particles.hpp @@ -125,6 +125,8 @@ class Particles // Simulation total domain. std::array global_mesh_ext; + std::array global_mesh_lo; + std::array global_mesh_hi; // Simulation sub domain (single MPI rank). std::array local_mesh_ext; @@ -201,6 +203,8 @@ class Particles std::array is_periodic; for ( int d = 0; d < dim; d++ ) { + global_mesh_lo[d] = global_mesh->lowCorner( d ); + global_mesh_hi[d] = global_mesh->highCorner( d ); global_mesh_ext[d] = global_mesh->extent( d ); is_periodic[d] = false; } @@ -516,6 +520,8 @@ class Particles // Simulation total domain. using base_type::global_mesh_ext; + using base_type::global_mesh_hi; + using base_type::global_mesh_lo; // Simulation sub domain (single MPI rank). using base_type::ghost_mesh_hi;