Compadre 1.5.5
Loading...
Searching...
No Matches
GMLS_Staggered_Manifold.cpp
Go to the documentation of this file.
1#include <iostream>
2#include <string>
3#include <vector>
4#include <map>
5#include <stdlib.h>
6#include <cstdio>
7#include <random>
8
9#include <Compadre_Config.h>
10#include <Compadre_GMLS.hpp>
13
14#include "GMLS_Manifold.hpp"
16
17#ifdef COMPADRE_USE_MPI
18#include <mpi.h>
19#endif
20
21#include <Kokkos_Timer.hpp>
22#include <Kokkos_Core.hpp>
23
24using namespace Compadre;
25
26//! [Parse Command Line Arguments]
27
28// called from command line
29int main (int argc, char* args[]) {
30
31// initializes MPI (if available) with command line arguments given
32#ifdef COMPADRE_USE_MPI
33MPI_Init(&argc, &args);
34#endif
35
36// initializes Kokkos with command line arguments given
37Kokkos::initialize(argc, args);
38
39// code block to reduce scope for all Kokkos View allocations
40// otherwise, Views may be deallocating when we call Kokkos::finalize() later
41{
42
43 CommandLineProcessor clp(argc, args);
44 auto order = clp.order;
45 auto dimension = clp.dimension;
46 auto number_target_coords = clp.number_target_coords;
47 auto constraint_name = clp.constraint_name;
48 auto solver_name = clp.solver_name;
49 auto problem_name = clp.problem_name;
50 int N_pts_on_sphere = (clp.number_source_coords>=0) ? clp.number_source_coords : 1000;
51
52 // minimum neighbors for unisolvency is the same as the size of the polynomial basis
53 // dimension has one subtracted because it is a D-1 manifold represented in D dimensions
54 const int min_neighbors = Compadre::GMLS::getNP(order, dimension-1);
55
56 //! [Parse Command Line Arguments]
57 Kokkos::Timer timer;
58 Kokkos::Profiling::pushRegion("Setup Point Data");
59 //! [Setting Up The Point Cloud]
60
61
62 // coordinates of source sites, bigger than needed then resized later
63 Kokkos::View<double**, Kokkos::DefaultExecutionSpace> source_coords_device("source coordinates",
64 1.25*N_pts_on_sphere, 3);
65 Kokkos::View<double**>::HostMirror source_coords = Kokkos::create_mirror_view(source_coords_device);
66
67 double r = 1.0;
68
69 // set number of source coordinates from what was calculated
70 int number_source_coords;
71
72 { // fill source coordinates from surface of a sphere with quasiuniform points
73 // algorithm described at https://www.cmu.edu/biolphys/deserno/pdf/sphere_equi.pdf
74 int N_count = 0;
75 double a = 4*PI*r*r/N_pts_on_sphere;
76 double d = std::sqrt(a);
77 int M_theta = std::round(PI/d);
78 double d_theta = PI/M_theta;
79 double d_phi = a/d_theta;
80 for (int i=0; i<M_theta; ++i) {
81 double theta = PI*(i + 0.5)/M_theta;
82 int M_phi = std::round(2*PI*std::sin(theta)/d_phi);
83 for (int j=0; j<M_phi; ++j) {
84 double phi = 2*PI*j/M_phi;
85 source_coords(N_count, 0) = theta;
86 source_coords(N_count, 1) = phi;
87 N_count++;
88 }
89 }
90
91 for (int i=0; i<N_count; ++i) {
92 double theta = source_coords(i,0);
93 double phi = source_coords(i,1);
94 source_coords(i,0) = r*std::sin(theta)*std::cos(phi);
95 source_coords(i,1) = r*std::sin(theta)*std::sin(phi);
96 source_coords(i,2) = r*cos(theta);
97 //printf("%f %f %f\n", source_coords(i,0), source_coords(i,1), source_coords(i,2));
98 }
99 number_source_coords = N_count;
100 }
101
102 // resize source_coords to the size actually needed
103 Kokkos::resize(source_coords, number_source_coords, 3);
104 Kokkos::resize(source_coords_device, number_source_coords, 3);
105
106 // coordinates of target sites
107 Kokkos::View<double**, Kokkos::DefaultExecutionSpace> target_coords_device("target coordinates",
108 number_target_coords, 3);
109 Kokkos::View<double**>::HostMirror target_coords = Kokkos::create_mirror_view(target_coords_device);
110
111 // seed random number generator
112 std::mt19937 rng(50);
113
114 // generate random integers in [0...number_source_coords-1] (used to pick target sites)
115 std::uniform_int_distribution<int> gen_num_neighbors(0, number_source_coords-1); // uniform, unbiased
116
117 // fill target sites with random selections from source sites
118 for (int i=0; i<number_target_coords; ++i) {
119 const int source_site_to_copy = gen_num_neighbors(rng);
120 for (int j=0; j<3; ++j) {
121 target_coords(i,j) = source_coords(source_site_to_copy,j);
122 }
123 }
124
125
126 //! [Setting Up The Point Cloud]
127
128 Kokkos::Profiling::popRegion();
129 Kokkos::fence();
130 Kokkos::Profiling::pushRegion("Creating Data");
131
132 //! [Creating The Data]
133
134
135 // source coordinates need copied to device before using to construct sampling data
136 Kokkos::deep_copy(source_coords_device, source_coords);
137 Kokkos::deep_copy(target_coords_device, target_coords);
138
139 // ensure that source coordinates are sent to device before evaluating sampling data based on them
140 Kokkos::fence();
141
142
143 // need Kokkos View storing true solution (for samples)
144 Kokkos::View<double*, Kokkos::DefaultExecutionSpace> sampling_data_device("samples of true solution",
145 source_coords_device.extent(0));
146
147 // need Kokkos View storing true vector solution (for samples)
148 Kokkos::View<double**, Kokkos::DefaultExecutionSpace> sampling_vector_data_device("samples of vector true solution",
149 source_coords_device.extent(0), 3);
150
151 Kokkos::parallel_for("Sampling Manufactured Solutions", Kokkos::RangePolicy<Kokkos::DefaultExecutionSpace>
152 (0,source_coords.extent(0)), KOKKOS_LAMBDA(const int i) {
153
154 // coordinates of source site i
155 double xval = source_coords_device(i,0);
156 double yval = (dimension>1) ? source_coords_device(i,1) : 0;
157 double zval = (dimension>2) ? source_coords_device(i,2) : 0;
158
159 // data for targets with scalar input
160 sampling_data_device(i) = sphere_harmonic54(xval, yval, zval);
161 //printf("%f\n", sampling_data_device(i));
162
163 for (int j=0; j<3; ++j) {
164 double gradient[3] = {0,0,0};
165 gradient_sphereHarmonic54_ambient(gradient, xval, yval, zval);
166 sampling_vector_data_device(i,j) = gradient[j];
167 }
168 //printf("%f %f %f\n", sampling_vector_data_device(i,0), sampling_vector_data_device(i,1), sampling_vector_data_device(i,2));
169 });
170
171
172 //! [Creating The Data]
173
174 Kokkos::Profiling::popRegion();
175 Kokkos::Profiling::pushRegion("Neighbor Search");
176
177 //! [Performing Neighbor Search]
178
179
180 // Point cloud construction for neighbor search
181 // CreatePointCloudSearch constructs an object of type PointCloudSearch, but deduces the templates for you
182 auto point_cloud_search(CreatePointCloudSearch(source_coords, dimension));
183
184 // each row is a neighbor list for a target site, with the first column of each row containing
185 // the number of neighbors for that rows corresponding target site
186 double epsilon_multiplier = 1.5;
187 int estimated_upper_bound_number_neighbors =
188 point_cloud_search.getEstimatedNumberNeighborsUpperBound(min_neighbors, dimension, epsilon_multiplier);
189
190 Kokkos::View<int**, Kokkos::DefaultExecutionSpace> neighbor_lists_device("neighbor lists",
191 number_target_coords, estimated_upper_bound_number_neighbors); // first column is # of neighbors
192 Kokkos::View<int**>::HostMirror neighbor_lists = Kokkos::create_mirror_view(neighbor_lists_device);
193
194 // each target site has a window size
195 Kokkos::View<double*, Kokkos::DefaultExecutionSpace> epsilon_device("h supports", number_target_coords);
196 Kokkos::View<double*>::HostMirror epsilon = Kokkos::create_mirror_view(epsilon_device);
197
198 // query the point cloud to generate the neighbor lists using a kdtree to produce the n nearest neighbor
199 // to each target site, adding (epsilon_multiplier-1)*100% to whatever the distance away the further neighbor used is from
200 // each target to the view for epsilon
201 point_cloud_search.generate2DNeighborListsFromKNNSearch(false /*not dry run*/, target_coords, neighbor_lists,
202 epsilon, min_neighbors, epsilon_multiplier);
203
204
205 //! [Performing Neighbor Search]
206
207 Kokkos::Profiling::popRegion();
208 Kokkos::fence(); // let call to build neighbor lists complete before copying back to device
209 timer.reset();
210
211 //! [Setting Up The GMLS Object]
212
213
214 // Copy data back to device (they were filled on the host)
215 // We could have filled Kokkos Views with memory space on the host
216 // and used these instead, and then the copying of data to the device
217 // would be performed in the GMLS class
218 Kokkos::deep_copy(neighbor_lists_device, neighbor_lists);
219 Kokkos::deep_copy(epsilon_device, epsilon);
220
221 // initialize an instance of the GMLS class
222 GMLS my_GMLS_vector_1(ReconstructionSpace::VectorTaylorPolynomial,
224 order, dimension,
225 solver_name.c_str(), problem_name.c_str(), constraint_name.c_str(),
226 order /*manifold order*/);
227
228 // pass in neighbor lists, source coordinates, target coordinates, and window sizes
229 //
230 // neighbor lists have the format:
231 // dimensions: (# number of target sites) X (# maximum number of neighbors for any given target + 1)
232 // the first column contains the number of neighbors for that rows corresponding target index
233 //
234 // source coordinates have the format:
235 // dimensions: (# number of source sites) X (dimension)
236 // entries in the neighbor lists (integers) correspond to rows of this 2D array
237 //
238 // target coordinates have the format:
239 // dimensions: (# number of target sites) X (dimension)
240 // # of target sites is same as # of rows of neighbor lists
241 //
242 my_GMLS_vector_1.setProblemData(neighbor_lists_device, source_coords_device, target_coords_device, epsilon_device);
243
244 // create a vector of target operations
245 std::vector<TargetOperation> lro_vector_1(1);
246 lro_vector_1[0] = DivergenceOfVectorPointEvaluation;
247
248 // and then pass them to the GMLS class
249 my_GMLS_vector_1.addTargets(lro_vector_1);
250
251 // sets the weighting kernel function from WeightingFunctionType for curvature
252 my_GMLS_vector_1.setCurvatureWeightingType(WeightingFunctionType::Power);
253
254 // power to use in the weighting kernel function for curvature coefficients
255 my_GMLS_vector_1.setCurvatureWeightingParameter(2);
256
257 // sets the weighting kernel function from WeightingFunctionType
258 my_GMLS_vector_1.setWeightingType(WeightingFunctionType::Power);
259
260 // power to use in that weighting kernel function
261 my_GMLS_vector_1.setWeightingParameter(2);
262
263 // setup quadrature for StaggeredEdgeIntegralSample
264 my_GMLS_vector_1.setOrderOfQuadraturePoints(2);
265 my_GMLS_vector_1.setDimensionOfQuadraturePoints(1);
266 my_GMLS_vector_1.setQuadratureType("LINE");
267
268 // generate the alphas that to be combined with data for each target operation requested in lro
269 my_GMLS_vector_1.generateAlphas();
270
271 // initialize another instance of the GMLS class
272 GMLS my_GMLS_vector_2(ReconstructionSpace::VectorTaylorPolynomial,
275 order, dimension,
276 solver_name.c_str(), problem_name.c_str(), constraint_name.c_str(),
277 order /*manifold order*/);
278
279 my_GMLS_vector_2.setProblemData(neighbor_lists_device, source_coords_device, target_coords_device, epsilon_device);
280 std::vector<TargetOperation> lro_vector_2(2);
282 lro_vector_2[1] = DivergenceOfVectorPointEvaluation;
283 //lro_vector_2[2] = GradientOfScalarPointEvaluation;
284 my_GMLS_vector_2.addTargets(lro_vector_2);
285 my_GMLS_vector_2.setCurvatureWeightingType(WeightingFunctionType::Power);
286 my_GMLS_vector_2.setCurvatureWeightingParameter(2);
287 my_GMLS_vector_2.setWeightingType(WeightingFunctionType::Power);
288 my_GMLS_vector_2.setWeightingParameter(2);
289 my_GMLS_vector_2.setOrderOfQuadraturePoints(2);
290 my_GMLS_vector_2.setDimensionOfQuadraturePoints(1);
291 my_GMLS_vector_2.setQuadratureType("LINE");
292 my_GMLS_vector_2.generateAlphas();
293
294 // initialize another instance of the GMLS class
295 GMLS my_GMLS_scalar(ReconstructionSpace::ScalarTaylorPolynomial,
297 order, dimension,
298 solver_name.c_str(), problem_name.c_str(), constraint_name.c_str(),
299 order /*manifold order*/);
300
301 my_GMLS_scalar.setProblemData(neighbor_lists_device, source_coords_device, target_coords_device, epsilon_device);
302
303 std::vector<TargetOperation> lro_scalar(1);
305 //lro_scalar[1] = GradientOfScalarPointEvaluation;
306 my_GMLS_scalar.addTargets(lro_scalar);
307 my_GMLS_scalar.setCurvatureWeightingType(WeightingFunctionType::Power);
308 my_GMLS_scalar.setCurvatureWeightingParameter(2);
309 my_GMLS_scalar.setWeightingType(WeightingFunctionType::Power);
310 my_GMLS_scalar.setWeightingParameter(2);
311 my_GMLS_scalar.generateAlphas();
312
313
314 //! [Setting Up The GMLS Object]
315
316 double instantiation_time = timer.seconds();
317 std::cout << "Took " << instantiation_time << "s to complete alphas generation." << std::endl;
318 Kokkos::fence(); // let generateAlphas finish up before using alphas
319 Kokkos::Profiling::pushRegion("Apply Alphas to Data");
320
321 //! [Apply GMLS Alphas To Data]
322
323
324 // it is important to note that if you expect to use the data as a 1D view, then you should use double*
325 // however, if you know that the target operation will result in a 2D view (vector or matrix output),
326 // then you should template with double** as this is something that can not be infered from the input data
327 // or the target operator at compile time. Additionally, a template argument is required indicating either
328 // Kokkos::HostSpace or Kokkos::DefaultExecutionSpace::memory_space()
329
330 // The Evaluator class takes care of handling input data views as well as the output data views.
331 // It uses information from the GMLS class to determine how many components are in the input
332 // as well as output for any choice of target functionals and then performs the contactions
333 // on the data using the alpha coefficients generated by the GMLS class, all on the device.
334 Evaluator vector_1_gmls_evaluator(&my_GMLS_vector_1);
335 Evaluator vector_2_gmls_evaluator(&my_GMLS_vector_2);
336 Evaluator scalar_gmls_evaluator(&my_GMLS_scalar);
337
338
339 //auto output_gradient_vectorbasis =
340 // vector_2_gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double**, Kokkos::HostSpace>
341 // (sampling_data_device, GradientOfScalarPointEvaluation, StaggeredEdgeAnalyticGradientIntegralSample);
342
343 //auto output_gradient_scalarbasis =
344 // scalar_gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double**, Kokkos::HostSpace>
345 // (sampling_data_device, GradientOfScalarPointEvaluation, StaggeredEdgeAnalyticGradientIntegralSample);
346
347 auto output_divergence_vectorsamples =
348 vector_1_gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double*, Kokkos::HostSpace>
350
351 auto output_divergence_scalarsamples =
352 vector_2_gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double*, Kokkos::HostSpace>
354
355 auto output_laplacian_vectorbasis =
356 vector_2_gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double*, Kokkos::HostSpace>
358
359 auto output_laplacian_scalarbasis =
360 scalar_gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double*, Kokkos::HostSpace>
362
363
364 //! [Apply GMLS Alphas To Data]
365
366 Kokkos::fence(); // let application of alphas to data finish before using results
367 Kokkos::Profiling::popRegion();
368 // times the Comparison in Kokkos
369 Kokkos::Profiling::pushRegion("Comparison");
370
371 //! [Check That Solutions Are Correct]
372
373 double laplacian_vectorbasis_error = 0;
374 double laplacian_vectorbasis_norm = 0;
375
376 double laplacian_scalarbasis_error = 0;
377 double laplacian_scalarbasis_norm = 0;
378
379 double gradient_vectorbasis_ambient_error = 0;
380 double gradient_vectorbasis_ambient_norm = 0;
381
382 double gradient_scalarbasis_ambient_error = 0;
383 double gradient_scalarbasis_ambient_norm = 0;
384
385 double divergence_vectorsamples_ambient_error = 0;
386 double divergence_vectorsamples_ambient_norm = 0;
387
388 double divergence_scalarsamples_ambient_error = 0;
389 double divergence_scalarsamples_ambient_norm = 0;
390
391 // loop through the target sites
392 for (int i=0; i<number_target_coords; i++) {
393
394 // target site i's coordinate
395 double xval = target_coords(i,0);
396 double yval = (dimension>1) ? target_coords(i,1) : 0;
397 double zval = (dimension>2) ? target_coords(i,2) : 0;
398
399 // evaluation of various exact solutions
400 double actual_Laplacian = laplace_beltrami_sphere_harmonic54(xval, yval, zval);
401 double actual_Gradient_ambient[3] = {0,0,0}; // initialized for 3, but only filled up to dimension
402 gradient_sphereHarmonic54_ambient(actual_Gradient_ambient, xval, yval, zval);
403
404 laplacian_vectorbasis_error += (output_laplacian_vectorbasis(i) - actual_Laplacian)*(output_laplacian_vectorbasis(i) - actual_Laplacian);
405 laplacian_vectorbasis_norm += actual_Laplacian*actual_Laplacian;
406
407 //printf("Error of %f, %f vs %f\n", (output_laplacian_scalarbasis(i) - actual_Laplacian), output_laplacian_scalarbasis(i), actual_Laplacian);
408 laplacian_scalarbasis_error += (output_laplacian_scalarbasis(i) - actual_Laplacian)*(output_laplacian_scalarbasis(i) - actual_Laplacian);
409 laplacian_scalarbasis_norm += actual_Laplacian*actual_Laplacian;
410
411 //for (int j=0; j<dimension; ++j) {
412 // //printf("VectorBasis Error of %f, %f vs %f\n", (output_gradient_vectorbasis(i,j) - actual_Gradient_ambient[j]), output_gradient_vectorbasis(i,j), actual_Gradient_ambient[j]);
413 // gradient_vectorbasis_ambient_error += (output_gradient_vectorbasis(i,j) - actual_Gradient_ambient[j])*(output_gradient_vectorbasis(i,j) - actual_Gradient_ambient[j]);
414 // gradient_vectorbasis_ambient_norm += actual_Gradient_ambient[j]*actual_Gradient_ambient[j];
415 //}
416
417 //for (int j=0; j<dimension; ++j) {
418 // //printf("ScalarBasis Error of %f, %f vs %f\n", (output_gradient_scalarbasis(i,j) - actual_Gradient_ambient[j]), output_gradient_scalarbasis(i,j), actual_Gradient_ambient[j]);
419 // gradient_scalarbasis_ambient_error += (output_gradient_scalarbasis(i,j) - actual_Gradient_ambient[j])*(output_gradient_scalarbasis(i,j) - actual_Gradient_ambient[j]);
420 // gradient_scalarbasis_ambient_norm += actual_Gradient_ambient[j]*actual_Gradient_ambient[j];
421 //}
422
423 //printf("Error of %f, %f vs %f\n", (output_divergence(i) - actual_Laplacian), output_divergence(i), actual_Laplacian);
424 divergence_vectorsamples_ambient_error += (output_divergence_vectorsamples(i) - actual_Laplacian)*(output_divergence_vectorsamples(i) - actual_Laplacian);
425 divergence_vectorsamples_ambient_norm += actual_Laplacian*actual_Laplacian;
426
427 divergence_scalarsamples_ambient_error += (output_divergence_scalarsamples(i) - actual_Laplacian)*(output_divergence_scalarsamples(i) - actual_Laplacian);
428 divergence_scalarsamples_ambient_norm += actual_Laplacian*actual_Laplacian;
429
430 }
431
432 laplacian_vectorbasis_error /= number_target_coords;
433 laplacian_vectorbasis_error = std::sqrt(laplacian_vectorbasis_error);
434 laplacian_vectorbasis_norm /= number_target_coords;
435 laplacian_vectorbasis_norm = std::sqrt(laplacian_vectorbasis_norm);
436
437 laplacian_scalarbasis_error /= number_target_coords;
438 laplacian_scalarbasis_error = std::sqrt(laplacian_scalarbasis_error);
439 laplacian_scalarbasis_norm /= number_target_coords;
440 laplacian_scalarbasis_norm = std::sqrt(laplacian_scalarbasis_norm);
441
442 gradient_vectorbasis_ambient_error /= number_target_coords;
443 gradient_vectorbasis_ambient_error = std::sqrt(gradient_vectorbasis_ambient_error);
444 gradient_vectorbasis_ambient_norm /= number_target_coords;
445 gradient_vectorbasis_ambient_norm = std::sqrt(gradient_vectorbasis_ambient_norm);
446
447 gradient_scalarbasis_ambient_error /= number_target_coords;
448 gradient_scalarbasis_ambient_error = std::sqrt(gradient_scalarbasis_ambient_error);
449 gradient_scalarbasis_ambient_norm /= number_target_coords;
450 gradient_scalarbasis_ambient_norm = std::sqrt(gradient_scalarbasis_ambient_norm);
451
452 divergence_vectorsamples_ambient_error /= number_target_coords;
453 divergence_vectorsamples_ambient_error = std::sqrt(divergence_vectorsamples_ambient_error);
454 divergence_vectorsamples_ambient_norm /= number_target_coords;
455 divergence_vectorsamples_ambient_norm = std::sqrt(divergence_vectorsamples_ambient_norm);
456
457 divergence_scalarsamples_ambient_error /= number_target_coords;
458 divergence_scalarsamples_ambient_error = std::sqrt(divergence_scalarsamples_ambient_error);
459 divergence_scalarsamples_ambient_norm /= number_target_coords;
460 divergence_scalarsamples_ambient_norm = std::sqrt(divergence_scalarsamples_ambient_norm);
461
462 printf("Staggered Laplace-Beltrami (VectorBasis) Error: %g\n", laplacian_vectorbasis_error / laplacian_vectorbasis_norm);
463 printf("Staggered Laplace-Beltrami (ScalarBasis) Error: %g\n", laplacian_scalarbasis_error / laplacian_scalarbasis_norm);
464 printf("Surface Staggered Gradient (VectorBasis) Error: %g\n", gradient_vectorbasis_ambient_error / gradient_vectorbasis_ambient_norm);
465 printf("Surface Staggered Gradient (ScalarBasis) Error: %g\n", gradient_scalarbasis_ambient_error / gradient_scalarbasis_ambient_norm);
466 printf("Surface Staggered Divergence (VectorSamples) Error: %g\n", divergence_vectorsamples_ambient_error / divergence_vectorsamples_ambient_norm);
467 printf("Surface Staggered Divergence (ScalarSamples) Error: %g\n", divergence_scalarsamples_ambient_error / divergence_scalarsamples_ambient_norm);
468 //! [Check That Solutions Are Correct]
469 // popRegion hidden from tutorial
470 // stop timing comparison loop
471 Kokkos::Profiling::popRegion();
472 //! [Finalize Program]
473
474
475} // end of code block to reduce scope, causing Kokkos View de-allocations
476// otherwise, Views may be deallocating when we call Kokkos::finalize() later
477
478// finalize Kokkos and MPI (if available)
479Kokkos::finalize();
480#ifdef COMPADRE_USE_MPI
481MPI_Finalize();
482#endif
483
484return 0;
485
486} // main
487
488
489//! [Finalize Program]
KOKKOS_INLINE_FUNCTION double sphere_harmonic54(double x, double y, double z)
#define PI
KOKKOS_INLINE_FUNCTION void gradient_sphereHarmonic54_ambient(double *gradient, double x, double y, double z)
KOKKOS_INLINE_FUNCTION double laplace_beltrami_sphere_harmonic54(double x, double y, double z)
int main(int argc, char *args[])
[Parse Command Line Arguments]
Lightweight Evaluator Helper This class is a lightweight wrapper for extracting and applying all rele...
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)
Generalized Moving Least Squares (GMLS)
void addTargets(TargetOperation lro)
Adds a target to the vector of target functional to be applied to the reconstruction.
void setDimensionOfQuadraturePoints(int dim)
Dimensions of quadrature points to use.
void setCurvatureWeightingParameter(int wp, int index=0)
Parameter for weighting kernel for curvature index = 0 sets p paramater for weighting kernel index = ...
void setWeightingParameter(int wp, int index=0)
Parameter for weighting kernel for GMLS problem index = 0 sets p paramater for weighting kernel index...
void generateAlphas(const int number_of_batches=1, const bool keep_coefficients=false, const bool clear_cache=true)
Meant to calculate target operations and apply the evaluations to the previously constructed polynomi...
void setOrderOfQuadraturePoints(int order)
Number quadrature points to use.
void setQuadratureType(std::string quadrature_type)
Type of quadrature points.
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)
void setWeightingType(const std::string &wt)
Type for weighting kernel for GMLS problem.
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....
void setCurvatureWeightingType(const std::string &wt)
Type for weighting kernel for curvature.
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...
constexpr SamplingFunctional StaggeredEdgeIntegralSample
Samples consist of the result of integrals of a vector dotted with the tangent along edges between ne...
@ ChainedStaggeredLaplacianOfScalarPointEvaluation
Point evaluation of the chained staggered Laplacian acting on VectorTaylorPolynomial basis + Staggere...
@ DivergenceOfVectorPointEvaluation
Point evaluation of the divergence of a vector (results in a scalar)
constexpr SamplingFunctional StaggeredEdgeAnalyticGradientIntegralSample
Analytical integral of a gradient source vector is just a difference of the scalar source at neighbor...