Compadre  1.6.4
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
GMLS_Device.cpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Compadre: COMpatible PArticle Discretization and REmap Toolkit
4 //
5 // Copyright 2018 NTESS and the Compadre contributors.
6 // SPDX-License-Identifier: BSD-2-Clause
7 // *****************************************************************************
8 // @HEADER
9 #include <iostream>
10 #include <string>
11 #include <vector>
12 #include <map>
13 #include <stdlib.h>
14 #include <cstdio>
15 #include <random>
16 
17 #include <Compadre_Config.h>
18 #include <Compadre_GMLS.hpp>
19 #include <Compadre_Evaluator.hpp>
21 
22 #include "GMLS_Tutorial.hpp"
23 #include "CommandLineProcessor.hpp"
24 
25 #ifdef COMPADRE_USE_MPI
26 #include <mpi.h>
27 #endif
28 
29 #include <Kokkos_Timer.hpp>
30 #include <Kokkos_Core.hpp>
31 
32 using namespace Compadre;
33 
34 //! [Parse Command Line Arguments]
35 
36 // called from command line
37 int main (int argc, char* args[]) {
38 
39 // initializes MPI (if available) with command line arguments given
40 #ifdef COMPADRE_USE_MPI
41 MPI_Init(&argc, &args);
42 #endif
43 
44 // initializes Kokkos with command line arguments given
45 Kokkos::initialize(argc, args);
46 
47 // becomes false if the computed solution not within the failure_threshold of the actual solution
48 bool all_passed = true;
49 
50 // code block to reduce scope for all Kokkos View allocations
51 // otherwise, Views may be deallocating when we call Kokkos finalize() later
52 {
53 
54  CommandLineProcessor clp(argc, args);
55  auto order = clp.order;
56  auto dimension = clp.dimension;
57  auto number_target_coords = clp.number_target_coords;
58  auto constraint_name = clp.constraint_name;
59  auto solver_name = clp.solver_name;
60  auto problem_name = clp.problem_name;
61  auto number_of_batches = clp.number_of_batches;
62  bool keep_coefficients = number_of_batches==1;
64 
65  // the functions we will be seeking to reconstruct are in the span of the basis
66  // of the reconstruction space we choose for GMLS, so the error should be very small
67  const double failure_tolerance = 1e-9;
68 
69  // Laplacian is a second order differential operator, which we expect to be slightly less accurate
70  const double laplacian_failure_tolerance = 1e-9;
71 
72  // minimum neighbors for unisolvency is the same as the size of the polynomial basis
73  const int min_neighbors = Compadre::GMLS::getNP(order, dimension);
74 
75  //! [Parse Command Line Arguments]
76  Kokkos::Timer timer;
77  Kokkos::Profiling::pushRegion("Setup Point Data");
78  //! [Setting Up The Point Cloud]
79 
80  // approximate spacing of source sites
81  double h_spacing = 0.05;
82  int n_neg1_to_1 = 2*(1/h_spacing) + 1; // always odd
83 
84  // number of source coordinate sites that will fill a box of [-1,1]x[-1,1]x[-1,1] with a spacing approximately h
85  const int number_source_coords = std::pow(n_neg1_to_1, dimension);
86 
87  // coordinates of source sites
88  Kokkos::View<double**, Kokkos::DefaultExecutionSpace> source_coords_device("source coordinates",
89  number_source_coords, 3);
90  Kokkos::View<double**>::host_mirror_type source_coords = Kokkos::create_mirror_view(source_coords_device);
91 
92  // coordinates of target sites
93  Kokkos::View<double**, Kokkos::DefaultExecutionSpace> target_coords_device ("target coordinates", number_target_coords, 3);
94  Kokkos::View<double**>::host_mirror_type target_coords = Kokkos::create_mirror_view(target_coords_device);
95 
96 
97  // fill source coordinates with a uniform grid
98  int source_index = 0;
99  double this_coord[3] = {0,0,0};
100  for (int i=-n_neg1_to_1/2; i<n_neg1_to_1/2+1; ++i) {
101  this_coord[0] = i*h_spacing;
102  for (int j=-n_neg1_to_1/2; j<n_neg1_to_1/2+1; ++j) {
103  this_coord[1] = j*h_spacing;
104  for (int k=-n_neg1_to_1/2; k<n_neg1_to_1/2+1; ++k) {
105  this_coord[2] = k*h_spacing;
106  if (dimension==3) {
107  source_coords(source_index,0) = this_coord[0];
108  source_coords(source_index,1) = this_coord[1];
109  source_coords(source_index,2) = this_coord[2];
110  source_index++;
111  }
112  }
113  if (dimension==2) {
114  source_coords(source_index,0) = this_coord[0];
115  source_coords(source_index,1) = this_coord[1];
116  source_coords(source_index,2) = 0;
117  source_index++;
118  }
119  }
120  if (dimension==1) {
121  source_coords(source_index,0) = this_coord[0];
122  source_coords(source_index,1) = 0;
123  source_coords(source_index,2) = 0;
124  source_index++;
125  }
126  }
127 
128  // fill target coords somewhere inside of [-0.5,0.5]x[-0.5,0.5]x[-0.5,0.5]
129  for(int i=0; i<number_target_coords; i++){
130 
131  // first, we get a uniformly random distributed direction
132  double rand_dir[3] = {0,0,0};
133 
134  for (int j=0; j<dimension; ++j) {
135  // rand_dir[j] is in [-0.5, 0.5]
136  rand_dir[j] = ((double)rand() / (double) RAND_MAX) - 0.5;
137  }
138 
139  // then we get a uniformly random radius
140  for (int j=0; j<dimension; ++j) {
141  target_coords(i,j) = rand_dir[j];
142  }
143 
144  }
145 
146 
147  //! [Setting Up The Point Cloud]
148 
149  Kokkos::Profiling::popRegion();
150  Kokkos::Profiling::pushRegion("Creating Data");
151 
152  //! [Creating The Data]
153 
154 
155  // source coordinates need copied to device before using to construct sampling data
156  Kokkos::deep_copy(source_coords_device, source_coords);
157 
158  // target coordinates copied next, because it is a convenient time to send them to device
159  Kokkos::deep_copy(target_coords_device, target_coords);
160 
161  // need Kokkos View storing true solution
162  Kokkos::View<double*, Kokkos::DefaultExecutionSpace> sampling_data_device("samples of true solution",
163  source_coords_device.extent(0));
164 
165  Kokkos::View<double**, Kokkos::DefaultExecutionSpace> gradient_sampling_data_device("samples of true gradient",
166  source_coords_device.extent(0), dimension);
167 
168  Kokkos::View<double**, Kokkos::DefaultExecutionSpace> divergence_sampling_data_device
169  ("samples of true solution for divergence test", source_coords_device.extent(0), dimension);
170 
171  Kokkos::parallel_for("Sampling Manufactured Solutions", Kokkos::RangePolicy<Kokkos::DefaultExecutionSpace>
172  (0,source_coords.extent(0)), KOKKOS_LAMBDA(const int i) {
173 
174  // coordinates of source site i
175  double xval = source_coords_device(i,0);
176  double yval = (dimension>1) ? source_coords_device(i,1) : 0;
177  double zval = (dimension>2) ? source_coords_device(i,2) : 0;
178 
179  // data for targets with scalar input
180  sampling_data_device(i) = trueSolution(xval, yval, zval, order, dimension);
181 
182  // data for targets with vector input (divergence)
183  double true_grad[3] = {0,0,0};
184  trueGradient(true_grad, xval, yval,zval, order, dimension);
185 
186  for (int j=0; j<dimension; ++j) {
187  gradient_sampling_data_device(i,j) = true_grad[j];
188 
189  // data for target with vector input (curl)
190  divergence_sampling_data_device(i,j) = divergenceTestSamples(xval, yval, zval, j, dimension);
191  }
192 
193  });
194 
195 
196  //! [Creating The Data]
197 
198  Kokkos::Profiling::popRegion();
199  Kokkos::Profiling::pushRegion("Neighbor Search");
200 
201  //! [Performing Neighbor Search]
202 
203 
204  // Point cloud construction for neighbor search
205  // CreatePointCloudSearch constructs an object of type PointCloudSearch, but deduces the templates for you
206  auto point_cloud_search(CreatePointCloudSearch(source_coords, dimension));
207 
208  double epsilon_multiplier = 1.4;
209 
210  // neighbor_lists_device will contain all neighbor lists (for each target site) in a compressed row format
211  // Initially, we do a dry-run to calculate neighborhood sizes before actually storing the result. This is
212  // why we can start with a neighbor_lists_device size of 0.
213  Kokkos::View<int*> neighbor_lists_device("neighbor lists",
214  0); // first column is # of neighbors
215  Kokkos::View<int*>::host_mirror_type neighbor_lists = Kokkos::create_mirror_view(neighbor_lists_device);
216  // number_of_neighbors_list must be the same size as the number of target sites so that it can be populated
217  // with the number of neighbors for each target site.
218  Kokkos::View<int*> number_of_neighbors_list_device("number of neighbor lists",
219  number_target_coords); // first column is # of neighbors
220  Kokkos::View<int*>::host_mirror_type number_of_neighbors_list = Kokkos::create_mirror_view(number_of_neighbors_list_device);
221 
222  // each target site has a window size
223  Kokkos::View<double*, Kokkos::DefaultExecutionSpace> epsilon_device("h supports", number_target_coords);
224  Kokkos::View<double*>::host_mirror_type epsilon = Kokkos::create_mirror_view(epsilon_device);
225 
226  // query the point cloud to generate the neighbor lists using a kdtree to produce the n nearest neighbor
227  // to each target site, adding (epsilon_multiplier-1)*100% to whatever the distance away the further neighbor used is from
228  // each target to the view for epsilon
229  //
230  // This dry run populates number_of_neighbors_list with neighborhood sizes
231  size_t storage_size = point_cloud_search.generateCRNeighborListsFromKNNSearch(true /*dry run*/, target_coords, neighbor_lists,
232  number_of_neighbors_list, epsilon, min_neighbors, epsilon_multiplier);
233 
234  // resize neighbor_lists_device so as to be large enough to contain all neighborhoods
235  Kokkos::resize(neighbor_lists_device, storage_size);
236  neighbor_lists = Kokkos::create_mirror_view(neighbor_lists_device);
237 
238  // query the point cloud a second time, but this time storing results into neighbor_lists
239  point_cloud_search.generateCRNeighborListsFromKNNSearch(false /*not dry run*/, target_coords, neighbor_lists,
240  number_of_neighbors_list, epsilon, min_neighbors, epsilon_multiplier);
241 
242  //! [Performing Neighbor Search]
243 
244  Kokkos::Profiling::popRegion();
245  Kokkos::fence(); // let call to build neighbor lists complete before copying back to device
246  timer.reset();
247 
248  //! [Setting Up The GMLS Object]
249 
250 
251  // Copy data back to device (they were filled on the host)
252  // We could have filled Kokkos Views with memory space on the host
253  // and used these instead, and then the copying of data to the device
254  // would be performed in the GMLS class
255  Kokkos::deep_copy(neighbor_lists_device, neighbor_lists);
256  Kokkos::deep_copy(number_of_neighbors_list_device, number_of_neighbors_list);
257  Kokkos::deep_copy(epsilon_device, epsilon);
258 
259  // initialize an instance of the GMLS class
261  order, dimension,
262  solver_name.c_str(), problem_name.c_str(), constraint_name.c_str(),
263  2 /*manifold order*/);
264 
265  // pass in neighbor lists, number of neighbor lists, source coordinates, target coordinates, and window sizes
266  //
267  // neighbor lists has a compressed row format and is a 1D view:
268  // dimensions: ? (determined by neighbor search, but total size of the sum of the number of neighbors over all target sites)
269  //
270  // number of neighbors list is a 1D view:
271  // dimensions: (# number of target sites)
272  // each entry contains the number of neighbors for a target site
273  //
274  // source coordinates have the format:
275  // dimensions: (# number of source sites) X (dimension)
276  // entries in the neighbor lists (integers) correspond to rows of this 2D array
277  //
278  // target coordinates have the format:
279  // dimensions: (# number of target sites) X (dimension)
280  // # of target sites is same as # of rows of neighbor lists
281  //
282  my_GMLS.setProblemData(neighbor_lists_device, number_of_neighbors_list_device, source_coords_device, target_coords_device, epsilon_device);
283 
284  // create a vector of target operations
285  std::vector<TargetOperation> lro(5);
286  lro[0] = ScalarPointEvaluation;
291 
292  // and then pass them to the GMLS class
293  my_GMLS.addTargets(lro);
294 
295  // sets the weighting kernel function from WeightingFunctionType
296  my_GMLS.setWeightingType(wt);
297 
298  // power to use in that weighting kernel function
299  my_GMLS.setWeightingParameter(2);
300 
301  // generate the alphas that to be combined with data for each target operation requested in lro
302  my_GMLS.generateAlphas(number_of_batches, keep_coefficients /* keep polynomial coefficients, only needed for a test later in this program */);
303 
304 
305  //! [Setting Up The GMLS Object]
306 
307  double instantiation_time = timer.seconds();
308  std::cout << "Took " << instantiation_time << "s to complete alphas generation." << std::endl;
309  Kokkos::fence(); // let generateAlphas finish up before using alphas
310  Kokkos::Profiling::pushRegion("Apply Alphas to Data");
311 
312  //! [Apply GMLS Alphas To Data]
313 
314  // it is important to note that if you expect to use the data as a 1D view, then you should use double*
315  // however, if you know that the target operation will result in a 2D view (vector or matrix output),
316  // then you should template with double** as this is something that can not be infered from the input data
317  // or the target operator at compile time. Additionally, a template argument is required indicating either
318  // Kokkos::HostSpace or Kokkos::DefaultExecutionSpace::memory_space()
319 
320  // The Evaluator class takes care of handling input data views as well as the output data views.
321  // It uses information from the GMLS class to determine how many components are in the input
322  // as well as output for any choice of target functionals and then performs the contactions
323  // on the data using the alpha coefficients generated by the GMLS class, all on the device.
324  Evaluator gmls_evaluator(&my_GMLS);
325 
326  auto output_value = gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double*, Kokkos::HostSpace>
327  (sampling_data_device, ScalarPointEvaluation);
328 
329  auto output_laplacian = gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double*, Kokkos::HostSpace>
330  (sampling_data_device, LaplacianOfScalarPointEvaluation);
331 
332  auto output_gradient = gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double**, Kokkos::HostSpace>
333  (sampling_data_device, GradientOfScalarPointEvaluation);
334 
335  auto output_divergence = gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double*, Kokkos::HostSpace>
336  (gradient_sampling_data_device, DivergenceOfVectorPointEvaluation, VectorPointSample);
337 
338  auto output_curl = gmls_evaluator.applyAlphasToDataAllComponentsAllTargetSites<double**, Kokkos::HostSpace>
339  (divergence_sampling_data_device, CurlOfVectorPointEvaluation);
340 
341  // retrieves polynomial coefficients instead of remapped field
342  decltype(output_curl) scalar_coefficients;
343  if (number_of_batches==1)
344  scalar_coefficients =
345  gmls_evaluator.applyFullPolynomialCoefficientsBasisToDataAllComponents<double**, Kokkos::HostSpace>
346  (sampling_data_device);
347 
348  //! [Apply GMLS Alphas To Data]
349 
350  Kokkos::fence(); // let application of alphas to data finish before using results
351  Kokkos::Profiling::popRegion();
352  // times the Comparison in Kokkos
353  Kokkos::Profiling::pushRegion("Comparison");
354 
355  //! [Check That Solutions Are Correct]
356 
357 
358  // loop through the target sites
359  for (int i=0; i<number_target_coords; i++) {
360 
361  // load value from output
362  double GMLS_value = output_value(i);
363 
364  // load laplacian from output
365  double GMLS_Laplacian = output_laplacian(i);
366 
367  // load partial x from gradient
368  // this is a test that the scalar_coefficients 2d array returned hold valid entries
369  // scalar_coefficients(i,1)*1./epsilon(i) is equivalent to the target operation acting
370  // on the polynomials applied to the polynomial coefficients
371  double GMLS_GradX = (number_of_batches==1) ? scalar_coefficients(i,1)*1./epsilon(i) : output_gradient(i,0);
372 
373  // load partial y from gradient
374  double GMLS_GradY = (dimension>1) ? output_gradient(i,1) : 0;
375 
376  // load partial z from gradient
377  double GMLS_GradZ = (dimension>2) ? output_gradient(i,2) : 0;
378 
379  // load divergence from output
380  double GMLS_Divergence = output_divergence(i);
381 
382  // load curl from output
383  double GMLS_CurlX = (dimension>1) ? output_curl(i,0) : 0;
384  double GMLS_CurlY = (dimension>1) ? output_curl(i,1) : 0;
385  double GMLS_CurlZ = (dimension>2) ? output_curl(i,2) : 0;
386 
387 
388  // target site i's coordinate
389  double xval = target_coords(i,0);
390  double yval = (dimension>1) ? target_coords(i,1) : 0;
391  double zval = (dimension>2) ? target_coords(i,2) : 0;
392 
393  // evaluation of various exact solutions
394  double actual_value = trueSolution(xval, yval, zval, order, dimension);
395  double actual_Laplacian = trueLaplacian(xval, yval, zval, order, dimension);
396 
397  double actual_Gradient[3] = {0,0,0}; // initialized for 3, but only filled up to dimension
398  trueGradient(actual_Gradient, xval, yval, zval, order, dimension);
399 
400  double actual_Divergence;
401  actual_Divergence = trueLaplacian(xval, yval, zval, order, dimension);
402 
403  double actual_Curl[3] = {0,0,0}; // initialized for 3, but only filled up to dimension
404  // (and not at all for dimimension = 1)
405  if (dimension>1) {
406  actual_Curl[0] = curlTestSolution(xval, yval, zval, 0, dimension);
407  actual_Curl[1] = curlTestSolution(xval, yval, zval, 1, dimension);
408  if (dimension>2) {
409  actual_Curl[2] = curlTestSolution(xval, yval, zval, 2, dimension);
410  }
411  }
412 
413  // check actual function value
414  if(GMLS_value!=GMLS_value || std::abs(actual_value - GMLS_value) > failure_tolerance) {
415  all_passed = false;
416  std::cout << i << " Failed Actual by: " << std::abs(actual_value - GMLS_value) << std::endl;
417  }
418 
419  // check Laplacian
420  if(std::abs(actual_Laplacian - GMLS_Laplacian) > laplacian_failure_tolerance) {
421  all_passed = false;
422  std::cout << i <<" Failed Laplacian by: " << std::abs(actual_Laplacian - GMLS_Laplacian) << std::endl;
423  }
424 
425  // check gradient
426  if(std::abs(actual_Gradient[0] - GMLS_GradX) > failure_tolerance) {
427  all_passed = false;
428  std::cout << i << " Failed GradX by: " << std::abs(actual_Gradient[0] - GMLS_GradX) << std::endl;
429  if (dimension>1) {
430  if(std::abs(actual_Gradient[1] - GMLS_GradY) > failure_tolerance) {
431  all_passed = false;
432  std::cout << i << " Failed GradY by: " << std::abs(actual_Gradient[1] - GMLS_GradY) << std::endl;
433  }
434  }
435  if (dimension>2) {
436  if(std::abs(actual_Gradient[2] - GMLS_GradZ) > failure_tolerance) {
437  all_passed = false;
438  std::cout << i << " Failed GradZ by: " << std::abs(actual_Gradient[2] - GMLS_GradZ) << std::endl;
439  }
440  }
441  }
442 
443  // check divergence
444  if(std::abs(actual_Divergence - GMLS_Divergence) > failure_tolerance) {
445  all_passed = false;
446  std::cout << i << " Failed Divergence by: " << std::abs(actual_Divergence - GMLS_Divergence) << std::endl;
447  }
448 
449  // check curl
450  if (order > 2) { // reconstructed solution not in basis unless order greater than 2 used
451  double tmp_diff = 0;
452  if (dimension>1)
453  tmp_diff += std::abs(actual_Curl[0] - GMLS_CurlX) + std::abs(actual_Curl[1] - GMLS_CurlY);
454  if (dimension>2)
455  tmp_diff += std::abs(actual_Curl[2] - GMLS_CurlZ);
456  if(std::abs(tmp_diff) > failure_tolerance) {
457  all_passed = false;
458  std::cout << i << " Failed Curl by: " << std::abs(tmp_diff) << std::endl;
459  }
460  }
461  }
462 
463 
464  //! [Check That Solutions Are Correct]
465  // popRegion hidden from tutorial
466  // stop timing comparison loop
467  Kokkos::Profiling::popRegion();
468  //! [Finalize Program]
469 
470 
471 } // end of code block to reduce scope, causing Kokkos View de-allocations
472 // otherwise, Views may be deallocating when we call Kokkos finalize() later
473 
474 // finalize Kokkos and MPI (if available)
475 Kokkos::finalize();
476 #ifdef COMPADRE_USE_MPI
477 MPI_Finalize();
478 #endif
479 
480 // output to user that test passed or failed
481 if(all_passed) {
482  fprintf(stdout, "Passed test \n");
483  return 0;
484 } else {
485  fprintf(stdout, "Failed test \n");
486  return -1;
487 }
488 
489 } // main
490 
491 
492 //! [Finalize Program]
Point evaluation of a scalar.
Scalar basis reused as many times as there are components in the vector resulting in a much cheaper p...
Lightweight Evaluator Helper This class is a lightweight wrapper for extracting and applying all rele...
Point evaluation of the gradient of a scalar.
WeightingFunctionType
Available weighting kernel function types.
Point evaluation of the laplacian of a scalar (could be on a manifold or not)
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...
int main(int argc, char **argv)
Point evaluation of the curl of a vector (results in a vector)
Kokkos::View< output_data_type, output_array_layout, output_memory_space > applyFullPolynomialCoefficientsBasisToDataAllComponents(view_type_input_data sampling_data, bool scalar_as_vector_if_needed=true) const
Generation of polynomial reconstruction coefficients by applying to data in GMLS (allocates memory fo...
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)
Point evaluation of the divergence of a vector (results in a scalar)
KOKKOS_INLINE_FUNCTION double trueLaplacian(double x, double y, double z, int order, int dimension)
Generalized Moving Least Squares (GMLS)
KOKKOS_INLINE_FUNCTION double curlTestSolution(double x, double y, double z, int component, int dimension)
KOKKOS_INLINE_FUNCTION double divergenceTestSamples(double x, double y, double z, int component, int dimension)
Compadre::WeightingFunctionType kernel_type
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)
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) ...
struct SubviewND< T, T2, enable_if_t<(T::rank< 2)> >{T _data_in;T2 _data_original_view;bool _scalar_as_vector_if_needed;SubviewND(T data_in, T2 data_original_view, bool scalar_as_vector_if_needed){_data_in=data_in;_data_original_view=data_original_view;_scalar_as_vector_if_needed=scalar_as_vector_if_needed;}auto get1DView(const int column_num) -> decltype(Kokkos::subview(_data_in, Kokkos::ALL))
Creates 1D subviews of data from a 1D view, generally constructed with CreateNDSliceOnDeviceView.
KOKKOS_INLINE_FUNCTION double trueSolution(double x, double y, double z, int order, int dimension)
constexpr SamplingFunctional VectorPointSample
Point evaluations of the entire vector source function.