Compadre  1.5.7
Compadre_BernsteinPolynomial.hpp
Go to the documentation of this file.
1 #ifndef _COMPADRE_GMLS_BERNSTEIN_BASIS_HPP_
2 #define _COMPADRE_GMLS_BERNSTEIN_BASIS_HPP_
3 
4 #include "Compadre_GMLS.hpp"
6 
7 namespace Compadre {
8 /*! \brief Definition of scalar Bernstein polynomial basis
9 
10  For order P, the sum of the basis is defined by:
11  \li in 1D) \f[\sum_{n=0}^{n=P} \frac{P!}{n!(P-n)!}*\left(\frac{x}{2h}+\frac{1}{2}\right)^n*\left(1-\frac{x}{2h}-\frac{1}{2}\right)^{P-n}\f]
12  \li in 2D) \f[\sum_{n_1=0}^{n_1=P}\sum_{n_2=0}^{n_2=P} \frac{P!}{n_1!(P-n_1)!}*\frac{P!}{n_2!(P-n_2)!}*\left(\frac{x}{2h}+\frac{1}{2}\right)^{n_1}*\left(1-\frac{x}{2h}-\frac{1}{2}\right)^{P-n_1}*\left(\frac{y}{2h}+\frac{1}{2}\right)^{n_2}*\left(1-\frac{y}{2h}-\frac{1}{2}\right)^{P-n_2}\f]
13  \li in 3D) \f[\sum_{n_1=0}^{n_1=P}\sum_{n_2=0}^{n_2=P}\sum_{n_3=0}^{n_3=P} \frac{P!}{n_1!(P-n_1)!}*\frac{P!}{n_2!(P-n_2)!}*\frac{P!}{n_3!(P-n_3)!}*\left(\frac{x}{2h}+\frac{1}{2}\right)^{n_1}*\left(1-\frac{x}{2h}-\frac{1}{2}\right)^{P-n_1}*\f]\f[\left(\frac{y}{2h}+\frac{1}{2}\right)^{n_2}*\left(1-\frac{y}{2h}-\frac{1}{2}\right)^{P-n_2}*\left(\frac{z}{2h}+\frac{1}{2}\right)^{n_3}*\left(1-\frac{z}{2h}-\frac{1}{2}\right)^{P-n_3}\f]
14 */
15 namespace BernsteinPolynomialBasis {
16 
17  /*! \brief Returns size of basis
18  \param degree [in] - highest degree of polynomial
19  \param dimension [in] - spatial dimension to evaluate
20  */
21  KOKKOS_INLINE_FUNCTION
22  int getSize(const int degree, const int dimension) {
23  return pown(ScalarTaylorPolynomialBasis::getSize(degree, int(1)), dimension);
24  }
25 
26  /*! \brief Evaluates the Bernstein polynomial basis
27  * delta[j] = weight_of_original_value * delta[j] + weight_of_new_value * (calculation of this function)
28  \param delta [in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _basis_multipler*the dimension of the polynomial basis.
29  \param workspace [in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
30  \param dimension [in] - spatial dimension to evaluate
31  \param max_degree [in] - highest degree of polynomial
32  \param h [in] - epsilon/window size
33  \param x [in] - x coordinate (already shifted by target)
34  \param y [in] - y coordinate (already shifted by target)
35  \param z [in] - z coordinate (already shifted by target)
36  \param starting_order [in] - polynomial order to start with (default=0)
37  \param weight_of_original_value [in] - weighting to assign to original value in delta (default=0, replace, 1-increment current value)
38  \param weight_of_new_value [in] - weighting to assign to new contribution (default=1, add to/replace)
39  */
40  KOKKOS_INLINE_FUNCTION
41  void evaluate(const member_type& teamMember, double* delta, double* workspace, const int dimension, const int max_degree, const int component, const double h, const double x, const double y, const double z, const int starting_order = 0, const double weight_of_original_value = 0.0, const double weight_of_new_value = 1.0) {
42  Kokkos::single(Kokkos::PerThread(teamMember), [&] () {
43  const double factorial[] = {1, 1, 2, 6, 24, 120, 720, 5040, 40320, 362880, 3628800, 39916800, 479001600, 6227020800, 87178291200};
44  if (dimension==3) {
45  scratch_vector_type x_over_h_to_i(workspace, max_degree+1);
46  scratch_vector_type y_over_h_to_i(workspace+1*(max_degree+1), max_degree+1);
47  scratch_vector_type z_over_h_to_i(workspace+2*(max_degree+1), max_degree+1);
48  scratch_vector_type one_minus_x_over_h_to_i(workspace+3*(max_degree+1), max_degree+1);
49  scratch_vector_type one_minus_y_over_h_to_i(workspace+4*(max_degree+1), max_degree+1);
50  scratch_vector_type one_minus_z_over_h_to_i(workspace+5*(max_degree+1), max_degree+1);
51  x_over_h_to_i[0] = 1;
52  y_over_h_to_i[0] = 1;
53  z_over_h_to_i[0] = 1;
54  one_minus_x_over_h_to_i[0] = 1;
55  one_minus_y_over_h_to_i[0] = 1;
56  one_minus_z_over_h_to_i[0] = 1;
57  for (int i=1; i<=max_degree; ++i) {
58  x_over_h_to_i[i] = x_over_h_to_i[i-1]*(0.5*x/h+0.5);
59  y_over_h_to_i[i] = y_over_h_to_i[i-1]*(0.5*y/h+0.5);
60  z_over_h_to_i[i] = z_over_h_to_i[i-1]*(0.5*z/h+0.5);
61  one_minus_x_over_h_to_i[i] = one_minus_x_over_h_to_i[i-1]*(1.0-(0.5*x/h+0.5));
62  one_minus_y_over_h_to_i[i] = one_minus_y_over_h_to_i[i-1]*(1.0-(0.5*y/h+0.5));
63  one_minus_z_over_h_to_i[i] = one_minus_z_over_h_to_i[i-1]*(1.0-(0.5*z/h+0.5));
64  }
65  // (in 3D) \sum_{n_1=0}^{n_1=P}\sum_{n_2=0}^{n_2=P}\sum_{n_3=0}^{n_3=P}
66  // \frac{P!}{n_1!(P-n_1)!}*\frac{P!}{n_2!(P-n_2)!}*\frac{P!}{n_3!(P-n_3)!}*
67  // \left(\frac{x}{2h}+\frac{1}{2}\right)^{n_1}*\left(1-\frac{x}{2h}-\frac{1}{2}\right)^{P-n_1}*
68  // \left(\frac{y}{2h}+\frac{1}{2}\right)^{n_2}*\left(1-\frac{y}{2h}-\frac{1}{2}\right)^{P-n_2}*
69  // \left(\frac{z}{2h}+\frac{1}{2}\right)^{n_3}*\left(1-\frac{z}{2h}-\frac{1}{2}\right)^{P-n_3}
70  int i=0;
71  for (int n1 = starting_order; n1 <= max_degree; n1++){
72  int degree_choose_n1 = factorial[max_degree]/(factorial[n1]*factorial[max_degree-n1]);
73  for (int n2 = starting_order; n2 <= max_degree; n2++){
74  int degree_choose_n2 = factorial[max_degree]/(factorial[n2]*factorial[max_degree-n2]);
75  for (int n3 = starting_order; n3 <= max_degree; n3++){
76  int degree_choose_n3 = factorial[max_degree]/(factorial[n3]*factorial[max_degree-n3]);
77  *(delta+i) = weight_of_original_value * *(delta+i) + weight_of_new_value * degree_choose_n1*degree_choose_n2*degree_choose_n3*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2]*z_over_h_to_i[n3]*one_minus_z_over_h_to_i[max_degree-n3];
78  i++;
79  }
80  }
81  }
82  } else if (dimension==2) {
83  scratch_vector_type x_over_h_to_i(workspace, max_degree+1);
84  scratch_vector_type y_over_h_to_i(workspace+1*(max_degree+1), max_degree+1);
85  scratch_vector_type one_minus_x_over_h_to_i(workspace+2*(max_degree+1), max_degree+1);
86  scratch_vector_type one_minus_y_over_h_to_i(workspace+3*(max_degree+1), max_degree+1);
87  x_over_h_to_i[0] = 1;
88  y_over_h_to_i[0] = 1;
89  one_minus_x_over_h_to_i[0] = 1;
90  one_minus_y_over_h_to_i[0] = 1;
91  for (int i=1; i<=max_degree; ++i) {
92  x_over_h_to_i[i] = x_over_h_to_i[i-1]*(0.5*x/h+0.5);
93  y_over_h_to_i[i] = y_over_h_to_i[i-1]*(0.5*y/h+0.5);
94  one_minus_x_over_h_to_i[i] = one_minus_x_over_h_to_i[i-1]*(1.0-(0.5*x/h+0.5));
95  one_minus_y_over_h_to_i[i] = one_minus_y_over_h_to_i[i-1]*(1.0-(0.5*y/h+0.5));
96  }
97  // (in 2D) \sum_{n_1=0}^{n_1=P}\sum_{n_2=0}^{n_2=P} \frac{P!}{n_1!(P-n_1)!}*
98  // \frac{P!}{n_2!(P-n_2)!}*\left(\frac{x}{2h}+\frac{1}{2}\right)^{n_1}*
99  // \left(1-\frac{x}{2h}-\frac{1}{2}\right)^{P-n_1}*\left(\frac{y}{2h}+\frac{1}{2}\right)^{n_2}*
100  // \left(1-\frac{y}{2h}-\frac{1}{2}\right)^{P-n_2}
101  int i = 0;
102  for (int n1 = starting_order; n1 <= max_degree; n1++){
103  int degree_choose_n1 = factorial[max_degree]/(factorial[n1]*factorial[max_degree-n1]);
104  for (int n2 = starting_order; n2 <= max_degree; n2++){
105  int degree_choose_n2 = factorial[max_degree]/(factorial[n2]*factorial[max_degree-n2]);
106  *(delta+i) = weight_of_original_value * *(delta+i) + weight_of_new_value * degree_choose_n1*degree_choose_n2*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2];
107  i++;
108  }
109  }
110  } else {
111  scratch_vector_type x_over_h_to_i(workspace, max_degree+1);
112  scratch_vector_type one_minus_x_over_h_to_i(workspace+1*(max_degree+1), max_degree+1);
113  x_over_h_to_i[0] = 1;
114  one_minus_x_over_h_to_i[0] = 1;
115  for (int i=1; i<=max_degree; ++i) {
116  x_over_h_to_i[i] = x_over_h_to_i[i-1]*(0.5*x/h+0.5);
117  one_minus_x_over_h_to_i[i] = one_minus_x_over_h_to_i[i-1]*(1.0-(0.5*x/h+0.5));
118  }
119  // (in 1D) \sum_{n=0}^{n=P} \frac{P!}{n!(P-n)!}*\left(\frac{x}{2h}+\frac{1}{2}\right)^n*
120  // \left(1-\frac{x}{2h}-\frac{1}{2}\right)^{P-n}
121  int i = 0;
122  for (int n=starting_order; n<=max_degree; ++n) {
123  int degree_choose_n = factorial[max_degree]/(factorial[n]*factorial[max_degree-n]);
124  *(delta+i) = weight_of_original_value * *(delta+i) + weight_of_new_value * degree_choose_n*x_over_h_to_i[n]*one_minus_x_over_h_to_i[max_degree-n];
125  i++;
126  }
127  }
128  });
129  }
130 
131  /*! \brief Evaluates the first partial derivatives of the Bernstein polynomial basis
132  * delta[j] = weight_of_original_value * delta[j] + weight_of_new_value * (calculation of this function)
133  \param delta [in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _basis_multipler*the dimension of the polynomial basis.
134  \param workspace [in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
135  \param dimension [in] - spatial dimension to evaluate
136  \param max_degree [in] - highest degree of polynomial
137  \param partial_direction [in] - direction in which to take partial derivative
138  \param h [in] - epsilon/window size
139  \param x [in] - x coordinate (already shifted by target)
140  \param y [in] - y coordinate (already shifted by target)
141  \param z [in] - z coordinate (already shifted by target)
142  \param starting_order [in] - polynomial order to start with (default=0)
143  \param weight_of_original_value [in] - weighting to assign to original value in delta (default=0, replace, 1-increment current value)
144  \param weight_of_new_value [in] - weighting to assign to new contribution (default=1, add to/replace)
145  */
146  KOKKOS_INLINE_FUNCTION
147  void evaluatePartialDerivative(const member_type& teamMember, double* delta, double* workspace, const int dimension, const int max_degree, const int component, const int partial_direction, const double h, const double x, const double y, const double z, const int starting_order = 0, const double weight_of_original_value = 0.0, const double weight_of_new_value = 1.0) {
148  Kokkos::single(Kokkos::PerThread(teamMember), [&] () {
149  const double factorial[] = {1, 1, 2, 6, 24, 120, 720, 5040, 40320, 362880, 3628800, 39916800, 479001600, 6227020800, 87178291200};
150  if (dimension==3) {
151  scratch_vector_type x_over_h_to_i(workspace, max_degree+1);
152  scratch_vector_type y_over_h_to_i(workspace+1*(max_degree+1), max_degree+1);
153  scratch_vector_type z_over_h_to_i(workspace+2*(max_degree+1), max_degree+1);
154  scratch_vector_type one_minus_x_over_h_to_i(workspace+3*(max_degree+1), max_degree+1);
155  scratch_vector_type one_minus_y_over_h_to_i(workspace+4*(max_degree+1), max_degree+1);
156  scratch_vector_type one_minus_z_over_h_to_i(workspace+5*(max_degree+1), max_degree+1);
157  x_over_h_to_i[0] = 1;
158  y_over_h_to_i[0] = 1;
159  z_over_h_to_i[0] = 1;
160  one_minus_x_over_h_to_i[0] = 1;
161  one_minus_y_over_h_to_i[0] = 1;
162  one_minus_z_over_h_to_i[0] = 1;
163  for (int i=1; i<=max_degree; ++i) {
164  x_over_h_to_i[i] = x_over_h_to_i[i-1]*(0.5*x/h+0.5);
165  y_over_h_to_i[i] = y_over_h_to_i[i-1]*(0.5*y/h+0.5);
166  z_over_h_to_i[i] = z_over_h_to_i[i-1]*(0.5*z/h+0.5);
167  one_minus_x_over_h_to_i[i] = one_minus_x_over_h_to_i[i-1]*(1.0-(0.5*x/h+0.5));
168  one_minus_y_over_h_to_i[i] = one_minus_y_over_h_to_i[i-1]*(1.0-(0.5*y/h+0.5));
169  one_minus_z_over_h_to_i[i] = one_minus_z_over_h_to_i[i-1]*(1.0-(0.5*z/h+0.5));
170  }
171  int i = 0;
172  for (int n1 = starting_order; n1 <= max_degree; n1++){
173  int degree_choose_n1 = factorial[max_degree]/(factorial[n1]*factorial[max_degree-n1]);
174  for (int n2 = starting_order; n2 <= max_degree; n2++){
175  int degree_choose_n2 = factorial[max_degree]/(factorial[n2]*factorial[max_degree-n2]);
176  for (int n3 = starting_order; n3 <= max_degree; n3++){
177  int degree_choose_n3 = factorial[max_degree]/(factorial[n3]*factorial[max_degree-n3]);
178  double new_contribution = 0.0;
179  if (partial_direction==0) {
180  if (n1>0) {
181  new_contribution += 0.5/h*n1*degree_choose_n1*degree_choose_n2*degree_choose_n3*x_over_h_to_i[n1-1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2]*z_over_h_to_i[n3]*one_minus_z_over_h_to_i[max_degree-n3];
182  }
183  if (max_degree-n1>0) {
184  new_contribution -= 0.5/h*(max_degree-n1)*degree_choose_n1*degree_choose_n2*degree_choose_n3*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1-1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2]*z_over_h_to_i[n3]*one_minus_z_over_h_to_i[max_degree-n3];
185  }
186  } else if (partial_direction==1) {
187  if (n2>0) {
188  new_contribution += 0.5/h*n2*degree_choose_n1*degree_choose_n2*degree_choose_n3*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2-1]*one_minus_y_over_h_to_i[max_degree-n2]*z_over_h_to_i[n3]*one_minus_z_over_h_to_i[max_degree-n3];
189  }
190  if (max_degree-n2>0) {
191  new_contribution -= 0.5/h*(max_degree-n2)*degree_choose_n1*degree_choose_n2*degree_choose_n3*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2-1]*z_over_h_to_i[n3]*one_minus_z_over_h_to_i[max_degree-n3];
192  }
193  } else {
194  if (n3>0) {
195  new_contribution += 0.5/h*n3*degree_choose_n1*degree_choose_n2*degree_choose_n3*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2]*z_over_h_to_i[n3-1]*one_minus_z_over_h_to_i[max_degree-n3];
196  }
197  if (max_degree-n3>0) {
198  new_contribution -= 0.5/h*(max_degree-n3)*degree_choose_n1*degree_choose_n2*degree_choose_n3*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2]*z_over_h_to_i[n3]*one_minus_z_over_h_to_i[max_degree-n3-1];
199  }
200  }
201  *(delta+i) = weight_of_original_value * *(delta+i) + weight_of_new_value * new_contribution;
202  i++;
203  }
204  }
205  }
206  } else if (dimension==2) {
207  scratch_vector_type x_over_h_to_i(workspace, max_degree+1);
208  scratch_vector_type y_over_h_to_i(workspace+1*(max_degree+1), max_degree+1);
209  scratch_vector_type one_minus_x_over_h_to_i(workspace+2*(max_degree+1), max_degree+1);
210  scratch_vector_type one_minus_y_over_h_to_i(workspace+3*(max_degree+1), max_degree+1);
211  x_over_h_to_i[0] = 1;
212  y_over_h_to_i[0] = 1;
213  one_minus_x_over_h_to_i[0] = 1;
214  one_minus_y_over_h_to_i[0] = 1;
215  for (int i=1; i<=max_degree; ++i) {
216  x_over_h_to_i[i] = x_over_h_to_i[i-1]*(0.5*x/h+0.5);
217  y_over_h_to_i[i] = y_over_h_to_i[i-1]*(0.5*y/h+0.5);
218  one_minus_x_over_h_to_i[i] = one_minus_x_over_h_to_i[i-1]*(1.0-(0.5*x/h+0.5));
219  one_minus_y_over_h_to_i[i] = one_minus_y_over_h_to_i[i-1]*(1.0-(0.5*y/h+0.5));
220  }
221  int i = 0;
222  for (int n1 = starting_order; n1 <= max_degree; n1++){
223  int degree_choose_n1 = factorial[max_degree]/(factorial[n1]*factorial[max_degree-n1]);
224  for (int n2 = starting_order; n2 <= max_degree; n2++){
225  int degree_choose_n2 = factorial[max_degree]/(factorial[n2]*factorial[max_degree-n2]);
226  double new_contribution = 0.0;
227  if (partial_direction==0) {
228  if (n1>0) {
229  new_contribution += 0.5/h*n1*degree_choose_n1*degree_choose_n2*x_over_h_to_i[n1-1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2];
230  }
231  if (max_degree-n1>0) {
232  new_contribution -= 0.5/h*(max_degree-n1)*degree_choose_n1*degree_choose_n2*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1-1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2];
233  }
234  } else {
235  if (n2>0) {
236  new_contribution += 0.5/h*n2*degree_choose_n1*degree_choose_n2*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2-1]*one_minus_y_over_h_to_i[max_degree-n2];
237  }
238  if (max_degree-n2>0) {
239  new_contribution -= 0.5/h*(max_degree-n2)*degree_choose_n1*degree_choose_n2*x_over_h_to_i[n1]*one_minus_x_over_h_to_i[max_degree-n1]*y_over_h_to_i[n2]*one_minus_y_over_h_to_i[max_degree-n2-1];
240  }
241  }
242  *(delta+i) = weight_of_original_value * *(delta+i) + weight_of_new_value * new_contribution;
243  i++;
244  }
245  }
246  } else {
247  scratch_vector_type x_over_h_to_i(workspace, max_degree+1);
248  scratch_vector_type one_minus_x_over_h_to_i(workspace+1*(max_degree+1), max_degree+1);
249  x_over_h_to_i[0] = 1;
250  one_minus_x_over_h_to_i[0] = 1;
251  for (int i=1; i<=max_degree; ++i) {
252  x_over_h_to_i[i] = x_over_h_to_i[i-1]*(0.5*x/h+0.5);
253  one_minus_x_over_h_to_i[i] = one_minus_x_over_h_to_i[i-1]*(1.0-(0.5*x/h+0.5));
254  }
255  int i = 0;
256  for (int n=starting_order; n<=max_degree; ++n) {
257  int degree_choose_n = factorial[max_degree]/(factorial[n]*factorial[max_degree-n]);
258  double new_contribution = 0.0;
259  if (n>0) {
260  //*(delta+i) = weight_of_original_value * *(delta+i) + weight_of_new_value * 0.5/h*n*degree_choose_n*x_over_h_to_i[n-1]*one_minus_x_over_h_to_i[max_degree-n];
261  new_contribution += 0.5/h*n*degree_choose_n*x_over_h_to_i[n-1]*one_minus_x_over_h_to_i[max_degree-n];
262  }
263  if (max_degree-n>0) {
264  //*(delta+i) = weight_of_original_value * *(delta+i) - weight_of_new_value * 0.5/h*(max_degree-n)*degree_choose_n*x_over_h_to_i[n]*one_minus_x_over_h_to_i[max_degree-n-1];
265  new_contribution -= 0.5/h*(max_degree-n)*degree_choose_n*x_over_h_to_i[n]*one_minus_x_over_h_to_i[max_degree-n-1];
266  }
267  *(delta+i) = weight_of_original_value * *(delta+i) + weight_of_new_value * new_contribution;
268  i++;
269  }
270  }
271  });
272  }
273 
274  /*! \brief Evaluates the second partial derivatives of the Bernstein polynomial basis
275  * delta[j] = weight_of_original_value * delta[j] + weight_of_new_value * (calculation of this function)
276  \param delta [in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _basis_multipler*the dimension of the polynomial basis.
277  \param workspace [in/out] - scratch space that is allocated so that each thread has its own copy. Must be at least as large as the _poly_order*the spatial dimension of the polynomial basis.
278  \param dimension [in] - spatial dimension to evaluate
279  \param max_degree [in] - highest degree of polynomial
280  \param partial_direction_1 [in] - direction in which to take first partial derivative
281  \param partial_direction_2 [in] - direction in which to take second partial derivative
282  \param h [in] - epsilon/window size
283  \param x [in] - x coordinate (already shifted by target)
284  \param y [in] - y coordinate (already shifted by target)
285  \param z [in] - z coordinate (already shifted by target)
286  \param starting_order [in] - polynomial order to start with (default=0)
287  \param weight_of_original_value [in] - weighting to assign to original value in delta (default=0, replace, 1-increment current value)
288  \param weight_of_new_value [in] - weighting to assign to new contribution (default=1, add to/replace)
289  */
290  KOKKOS_INLINE_FUNCTION
291  void evaluateSecondPartialDerivative(const member_type& teamMember, double* delta, double* workspace, const int dimension, const int max_degree, const int component, const int partial_direction_1, const int partial_direction_2, const double h, const double x, const double y, const double z, const int starting_order = 0, const double weight_of_original_value = 0.0, const double weight_of_new_value = 1.0) {
292  Kokkos::single(Kokkos::PerThread(teamMember), [&] () {
293  compadre_kernel_assert_release((false) && "Second partial derivatives not available for Bernstein basis.");
294  });
295  }
296 
297 } // BernsteinPolynomialBasis
298 } // Compadre
299 
300 #endif
Kokkos::View< double *, Kokkos::MemoryTraits< Kokkos::Unmanaged > > scratch_vector_type
#define compadre_kernel_assert_release(condition)
compadre_kernel_assert_release is similar to compadre_assert_release, but is a call on the device,...
team_policy::member_type member_type
KOKKOS_INLINE_FUNCTION void evaluatePartialDerivative(const member_type &teamMember, double *delta, double *workspace, const int dimension, const int max_degree, const int component, const int partial_direction, const double h, const double x, const double y, const double z, const int starting_order=0, const double weight_of_original_value=0.0, const double weight_of_new_value=1.0)
Evaluates the first partial derivatives of the Bernstein polynomial basis delta[j] = weight_of_origin...
KOKKOS_INLINE_FUNCTION void evaluateSecondPartialDerivative(const member_type &teamMember, double *delta, double *workspace, const int dimension, const int max_degree, const int component, const int partial_direction_1, const int partial_direction_2, const double h, const double x, const double y, const double z, const int starting_order=0, const double weight_of_original_value=0.0, const double weight_of_new_value=1.0)
Evaluates the second partial derivatives of the Bernstein polynomial basis delta[j] = weight_of_origi...
KOKKOS_INLINE_FUNCTION int getSize(const int degree, const int dimension)
Returns size of basis.
KOKKOS_INLINE_FUNCTION void evaluate(const member_type &teamMember, double *delta, double *workspace, const int dimension, const int max_degree, const int component, const double h, const double x, const double y, const double z, const int starting_order=0, const double weight_of_original_value=0.0, const double weight_of_new_value=1.0)
Evaluates the Bernstein polynomial basis delta[j] = weight_of_original_value * delta[j] + weight_of_n...
KOKKOS_INLINE_FUNCTION int getSize(const int degree, const int dimension)
Returns size of basis.
KOKKOS_INLINE_FUNCTION int pown(int n, unsigned p)
n^p (n^0 returns 1, regardless of n)