9 #include <Compadre_Config.h> 18 #ifdef COMPADRE_USE_MPI 22 #include <Kokkos_Timer.hpp> 23 #include <Kokkos_Core.hpp> 30 int main (
int argc,
char* args[]) {
33 #ifdef COMPADRE_USE_MPI 34 MPI_Init(&argc, &args);
41 bool all_passed =
true;
48 auto order = clp.
order;
58 const double failure_tolerance = 1e-9;
65 Kokkos::Profiling::pushRegion(
"Setup Point Data");
69 double h_spacing = 0.05;
70 int n_neg1_to_1 = 2*(1/h_spacing) + 1;
73 const int number_source_coords = std::pow(n_neg1_to_1, dimension);
76 Kokkos::View<double**, Kokkos::DefaultExecutionSpace> source_coords_device(
"source coordinates",
77 number_source_coords, 3);
78 Kokkos::View<double**>::HostMirror source_coords = Kokkos::create_mirror_view(source_coords_device);
81 Kokkos::View<double**, Kokkos::DefaultExecutionSpace> target_coords_device (
"target coordinates", number_target_coords, 3);
82 Kokkos::View<double**>::HostMirror target_coords = Kokkos::create_mirror_view(target_coords_device);
85 Kokkos::View<double***, Kokkos::DefaultExecutionSpace> tangent_bundles_device (
"tangent bundles", number_target_coords, dimension, dimension);
86 Kokkos::View<double***>::HostMirror tangent_bundles = Kokkos::create_mirror_view(tangent_bundles_device);
90 double this_coord[3] = {0,0,0};
91 for (
int i=-n_neg1_to_1/2; i<n_neg1_to_1/2+1; ++i) {
92 this_coord[0] = i*h_spacing;
93 for (
int j=-n_neg1_to_1/2; j<n_neg1_to_1/2+1; ++j) {
94 this_coord[1] = j*h_spacing;
95 for (
int k=-n_neg1_to_1/2; k<n_neg1_to_1/2+1; ++k) {
96 this_coord[2] = k*h_spacing;
98 source_coords(source_index,0) = this_coord[0];
99 source_coords(source_index,1) = this_coord[1];
100 source_coords(source_index,2) = this_coord[2];
105 source_coords(source_index,0) = this_coord[0];
106 source_coords(source_index,1) = this_coord[1];
107 source_coords(source_index,2) = 0;
112 source_coords(source_index,0) = this_coord[0];
113 source_coords(source_index,1) = 0;
114 source_coords(source_index,2) = 0;
120 for(
int i=0; i<number_target_coords; i++){
123 double rand_dir[3] = {0,0,0};
125 for (
int j=0; j<dimension; ++j) {
127 rand_dir[j] = ((double)rand() / (double) RAND_MAX) - 0.5;
131 for (
int j=0; j<dimension; ++j) {
132 target_coords(i,j) = rand_dir[j];
137 if (dimension == 2) {
138 tangent_bundles(i, 0, 0) = 0.0;
139 tangent_bundles(i, 0, 1) = 0.0;
140 tangent_bundles(i, 1, 0) = 1.0/(sqrt(2.0));
141 tangent_bundles(i, 1, 1) = 1.0/(sqrt(2.0));
142 }
else if (dimension == 3) {
143 tangent_bundles(i, 0, 0) = 0.0;
144 tangent_bundles(i, 0, 1) = 0.0;
145 tangent_bundles(i, 0, 2) = 0.0;
146 tangent_bundles(i, 1, 0) = 0.0;
147 tangent_bundles(i, 1, 1) = 0.0;
148 tangent_bundles(i, 1, 2) = 0.0;
149 tangent_bundles(i, 2, 0) = 1.0/(sqrt(3.0));
150 tangent_bundles(i, 2, 1) = 1.0/(sqrt(3.0));
151 tangent_bundles(i, 2, 2) = 1.0/(sqrt(3.0));
157 Kokkos::Profiling::popRegion();
158 Kokkos::Profiling::pushRegion(
"Creating Data");
164 Kokkos::deep_copy(source_coords_device, source_coords);
167 Kokkos::deep_copy(target_coords_device, target_coords);
170 Kokkos::deep_copy(tangent_bundles_device, tangent_bundles);
173 Kokkos::View<double*, Kokkos::DefaultExecutionSpace> sampling_data_device(
"samples of true solution",
174 source_coords_device.extent(0));
176 Kokkos::parallel_for(
"Sampling Manufactured Solutions", Kokkos::RangePolicy<Kokkos::DefaultExecutionSpace>
177 (0,source_coords.extent(0)), KOKKOS_LAMBDA(
const int i) {
180 double xval = source_coords_device(i,0);
181 double yval = (dimension>1) ? source_coords_device(i,1) : 0;
182 double zval = (dimension>2) ? source_coords_device(i,2) : 0;
185 sampling_data_device(i) =
trueSolution(xval, yval, zval, order, dimension);
190 Kokkos::Profiling::popRegion();
191 Kokkos::Profiling::pushRegion(
"Neighbor Search");
202 double epsilon_multiplier = 1.8;
203 int estimated_upper_bound_number_neighbors =
204 point_cloud_search.getEstimatedNumberNeighborsUpperBound(min_neighbors, dimension, epsilon_multiplier);
206 Kokkos::View<int**, Kokkos::DefaultExecutionSpace> neighbor_lists_device(
"neighbor lists",
207 number_target_coords, estimated_upper_bound_number_neighbors);
208 Kokkos::View<int**>::HostMirror neighbor_lists = Kokkos::create_mirror_view(neighbor_lists_device);
211 Kokkos::View<double*, Kokkos::DefaultExecutionSpace> epsilon_device(
"h supports", number_target_coords);
212 Kokkos::View<double*>::HostMirror epsilon = Kokkos::create_mirror_view(epsilon_device);
217 point_cloud_search.generate2DNeighborListsFromKNNSearch(
false , target_coords, neighbor_lists,
218 epsilon, min_neighbors, epsilon_multiplier);
222 Kokkos::Profiling::popRegion();
233 Kokkos::deep_copy(neighbor_lists_device, neighbor_lists);
234 Kokkos::deep_copy(epsilon_device, epsilon);
240 solver_name.c_str(), problem_name.c_str(), constraint_name.c_str(),
257 my_GMLS.
setProblemData(neighbor_lists_device, source_coords_device, target_coords_device, epsilon_device);
258 my_GMLS.setTangentBundle(tangent_bundles_device);
265 my_GMLS.addTargets(lro);
271 my_GMLS.setWeightingPower(2);
274 my_GMLS.generateAlphas(number_of_batches);
278 double instantiation_time = timer.seconds();
279 std::cout <<
"Took " << instantiation_time <<
"s to complete normal vectors generation." << std::endl;
281 Kokkos::Profiling::pushRegion(
"Apply Alphas to Data");
301 Kokkos::Profiling::popRegion();
303 Kokkos::Profiling::pushRegion(
"Comparison");
308 for (
int i=0; i<number_target_coords; i++) {
310 double xval = target_coords(i,0);
311 double yval = (dimension>1) ? target_coords(i,1) : 0;
312 double zval = (dimension>2) ? target_coords(i,2) : 0;
315 int num_neigh_i = neighbor_lists(i, 0);
316 double b_i = my_GMLS.getAlpha0TensorTo0Tensor(lro, i, num_neigh_i);
319 double GMLS_value = output_value(i);
322 double actual_Laplacian =
trueLaplacian(xval, yval, zval, order, dimension);
325 double actual_Gradient[3] = {0,0,0};
326 trueGradient(actual_Gradient, xval, yval, zval, order, dimension);
327 double g = (dimension == 3) ? (1.0/sqrt(3.0))*(actual_Gradient[0] + actual_Gradient[1] + actual_Gradient[2])
328 : (1.0/sqrt(2.0))*(actual_Gradient[0] + actual_Gradient[1]);
329 double adjusted_value = GMLS_value + b_i*g;
332 if(GMLS_value!=GMLS_value || std::abs(actual_Laplacian - adjusted_value) > failure_tolerance) {
334 std::cout << i <<
" Failed Actual by: " << std::abs(actual_Laplacian - adjusted_value) << std::endl;
341 Kokkos::Profiling::popRegion();
348 #ifdef COMPADRE_USE_MPI 354 fprintf(stdout,
"Passed test \n");
357 fprintf(stdout,
"Failed test \n");
Class handling Kokkos command line arguments and returning parameters.
Lightweight Evaluator Helper This class is a lightweight wrapper for extracting and applying all rele...
PointCloudSearch< view_type > CreatePointCloudSearch(view_type src_view, const local_index_type dimensions=-1, const local_index_type max_leaf=-1)
CreatePointCloudSearch allows for the construction of an object of type PointCloudSearch with templat...
Scalar polynomial basis centered at the target site and scaled by sum of basis powers e...
static KOKKOS_INLINE_FUNCTION int getNP(const int m, const int dimension=3, const ReconstructionSpace r_space=ReconstructionSpace::ScalarTaylorPolynomial)
Returns size of the basis for a given polynomial order and dimension General to dimension 1...
KOKKOS_INLINE_FUNCTION void trueGradient(double *ans, double x, double y, double z, int order, int dimension)
Kokkos::View< output_data_type, output_array_layout, output_memory_space > applyAlphasToDataAllComponentsAllTargetSites(view_type_input_data sampling_data, TargetOperation lro, const SamplingFunctional sro_in=PointSample, bool scalar_as_vector_if_needed=true, const int evaluation_site_local_index=0) const
Transformation of data under GMLS (allocates memory for output)
int main(int argc, char *args[])
[Parse Command Line Arguments]
Point evaluation of the laplacian of a scalar (could be on a manifold or not)
KOKKOS_INLINE_FUNCTION double trueLaplacian(double x, double y, double z, int order, int dimension)
Generalized Moving Least Squares (GMLS)
void setProblemData(view_type_1 neighbor_lists, view_type_2 source_coordinates, view_type_3 target_coordinates, view_type_4 epsilons)
Sets basic problem data (neighbor lists, source coordinates, and target coordinates) ...
constexpr SamplingFunctional PointSample
Available sampling functionals.
std::string constraint_name
KOKKOS_INLINE_FUNCTION double trueSolution(double x, double y, double z, int order, int dimension)
TargetOperation
Available target functionals.