Reference documentation for deal.II version Git 409ee4b167 2020-08-14 09:46:12 -0400
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
function_lib.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2019 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_function_lib_h
17 #define dealii_function_lib_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/function.h>
23 #include <deal.II/base/point.h>
25 #include <deal.II/base/table.h>
26 
27 #include <array>
28 
30 
39 namespace Functions
40 {
50  template <int dim>
51  class SquareFunction : public Function<dim>
52  {
53  public:
54  virtual double
55  value(const Point<dim> &p, const unsigned int component = 0) const override;
56  virtual void
57  vector_value(const Point<dim> &p, Vector<double> &values) const override;
58  virtual void
59  value_list(const std::vector<Point<dim>> &points,
60  std::vector<double> & values,
61  const unsigned int component = 0) const override;
62  virtual Tensor<1, dim>
63  gradient(const Point<dim> & p,
64  const unsigned int component = 0) const override;
65  virtual void
66  vector_gradient(const Point<dim> & p,
67  std::vector<Tensor<1, dim>> &gradient) const override;
68  virtual void
69  gradient_list(const std::vector<Point<dim>> &points,
70  std::vector<Tensor<1, dim>> & gradients,
71  const unsigned int component = 0) const override;
72  virtual double
73  laplacian(const Point<dim> & p,
74  const unsigned int component = 0) const override;
75  virtual void
76  laplacian_list(const std::vector<Point<dim>> &points,
77  std::vector<double> & values,
78  const unsigned int component = 0) const override;
79  };
80 
81 
82 
89  template <int dim>
90  class Q1WedgeFunction : public Function<dim>
91  {
92  public:
93  virtual double
94  value(const Point<dim> &p, const unsigned int component = 0) const override;
95 
96  virtual void
97  value_list(const std::vector<Point<dim>> &points,
98  std::vector<double> & values,
99  const unsigned int component = 0) const override;
100 
101  virtual void
102  vector_value_list(const std::vector<Point<dim>> &points,
103  std::vector<Vector<double>> & values) const override;
104 
105  virtual Tensor<1, dim>
106  gradient(const Point<dim> & p,
107  const unsigned int component = 0) const override;
108 
109  virtual void
110  gradient_list(const std::vector<Point<dim>> &points,
111  std::vector<Tensor<1, dim>> & gradients,
112  const unsigned int component = 0) const override;
113 
114  virtual void
116  const std::vector<Point<dim>> &,
117  std::vector<std::vector<Tensor<1, dim>>> &) const override;
118 
122  virtual double
123  laplacian(const Point<dim> & p,
124  const unsigned int component = 0) const override;
125 
129  virtual void
130  laplacian_list(const std::vector<Point<dim>> &points,
131  std::vector<double> & values,
132  const unsigned int component = 0) const override;
133  };
134 
135 
136 
151  template <int dim>
152  class PillowFunction : public Function<dim>
153  {
154  public:
159  PillowFunction(const double offset = 0.);
160 
164  virtual double
165  value(const Point<dim> &p, const unsigned int component = 0) const override;
166 
170  virtual void
171  value_list(const std::vector<Point<dim>> &points,
172  std::vector<double> & values,
173  const unsigned int component = 0) const override;
174 
178  virtual Tensor<1, dim>
179  gradient(const Point<dim> & p,
180  const unsigned int component = 0) const override;
181 
185  virtual void
186  gradient_list(const std::vector<Point<dim>> &points,
187  std::vector<Tensor<1, dim>> & gradients,
188  const unsigned int component = 0) const override;
189 
193  virtual double
194  laplacian(const Point<dim> & p,
195  const unsigned int component = 0) const override;
196 
200  virtual void
201  laplacian_list(const std::vector<Point<dim>> &points,
202  std::vector<double> & values,
203  const unsigned int component = 0) const override;
204 
205  private:
206  const double offset;
207  };
208 
209 
210 
218  template <int dim>
219  class CosineFunction : public Function<dim>
220  {
221  public:
226  CosineFunction(const unsigned int n_components = 1);
227 
228  virtual double
229  value(const Point<dim> &p, const unsigned int component = 0) const override;
230 
231  virtual void
232  value_list(const std::vector<Point<dim>> &points,
233  std::vector<double> & values,
234  const unsigned int component = 0) const override;
235 
236  virtual void
237  vector_value_list(const std::vector<Point<dim>> &points,
238  std::vector<Vector<double>> & values) const override;
239 
240  virtual Tensor<1, dim>
241  gradient(const Point<dim> & p,
242  const unsigned int component = 0) const override;
243 
244  virtual void
245  gradient_list(const std::vector<Point<dim>> &points,
246  std::vector<Tensor<1, dim>> & gradients,
247  const unsigned int component = 0) const override;
248 
249  virtual double
250  laplacian(const Point<dim> & p,
251  const unsigned int component = 0) const override;
252 
253  virtual void
254  laplacian_list(const std::vector<Point<dim>> &points,
255  std::vector<double> & values,
256  const unsigned int component = 0) const override;
257 
262  hessian(const Point<dim> & p,
263  const unsigned int component = 0) const override;
264 
268  virtual void
269  hessian_list(const std::vector<Point<dim>> & points,
270  std::vector<SymmetricTensor<2, dim>> &hessians,
271  const unsigned int component = 0) const override;
272  };
273 
274 
275 
286  template <int dim>
287  class CosineGradFunction : public Function<dim>
288  {
289  public:
294 
295  virtual double
296  value(const Point<dim> &p, const unsigned int component) const override;
297  virtual void
298  vector_value(const Point<dim> &p, Vector<double> &values) const override;
299  virtual void
300  value_list(const std::vector<Point<dim>> &points,
301  std::vector<double> & values,
302  const unsigned int component) const override;
303 
304  virtual void
305  vector_value_list(const std::vector<Point<dim>> &points,
306  std::vector<Vector<double>> & values) const override;
307 
308  virtual Tensor<1, dim>
309  gradient(const Point<dim> &p, const unsigned int component) const override;
310 
311  virtual void
312  gradient_list(const std::vector<Point<dim>> &points,
313  std::vector<Tensor<1, dim>> & gradients,
314  const unsigned int component) const override;
315 
316  virtual void
318  const std::vector<Point<dim>> & points,
319  std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;
320 
321  virtual double
322  laplacian(const Point<dim> &p, const unsigned int component) const override;
323  };
324 
325 
326 
332  template <int dim>
333  class ExpFunction : public Function<dim>
334  {
335  public:
339  virtual double
340  value(const Point<dim> &p, const unsigned int component = 0) const override;
341 
345  virtual void
346  value_list(const std::vector<Point<dim>> &points,
347  std::vector<double> & values,
348  const unsigned int component = 0) const override;
349 
353  virtual Tensor<1, dim>
354  gradient(const Point<dim> & p,
355  const unsigned int component = 0) const override;
356 
360  virtual void
361  gradient_list(const std::vector<Point<dim>> &points,
362  std::vector<Tensor<1, dim>> & gradients,
363  const unsigned int component = 0) const override;
364 
368  virtual double
369  laplacian(const Point<dim> & p,
370  const unsigned int component = 0) const override;
371 
375  virtual void
376  laplacian_list(const std::vector<Point<dim>> &points,
377  std::vector<double> & values,
378  const unsigned int component = 0) const override;
379  };
380 
381 
382 
392  class LSingularityFunction : public Function<2>
393  {
394  public:
395  virtual double
396  value(const Point<2> &p, const unsigned int component = 0) const override;
397 
398  virtual void
399  value_list(const std::vector<Point<2>> &points,
400  std::vector<double> & values,
401  const unsigned int component = 0) const override;
402 
403  virtual void
404  vector_value_list(const std::vector<Point<2>> &points,
405  std::vector<Vector<double>> &values) const override;
406 
407  virtual Tensor<1, 2>
408  gradient(const Point<2> & p,
409  const unsigned int component = 0) const override;
410 
411  virtual void
412  gradient_list(const std::vector<Point<2>> &points,
413  std::vector<Tensor<1, 2>> & gradients,
414  const unsigned int component = 0) const override;
415 
416  virtual void
418  const std::vector<Point<2>> &,
419  std::vector<std::vector<Tensor<1, 2>>> &) const override;
420 
421  virtual double
422  laplacian(const Point<2> & p,
423  const unsigned int component = 0) const override;
424 
425  virtual void
426  laplacian_list(const std::vector<Point<2>> &points,
427  std::vector<double> & values,
428  const unsigned int component = 0) const override;
429  };
430 
431 
432 
442  {
443  public:
448  virtual double
449  value(const Point<2> &p, const unsigned int component) const override;
450 
451  virtual void
452  value_list(const std::vector<Point<2>> &points,
453  std::vector<double> & values,
454  const unsigned int component) const override;
455 
456  virtual void
457  vector_value_list(const std::vector<Point<2>> &points,
458  std::vector<Vector<double>> &values) const override;
459 
460  virtual Tensor<1, 2>
461  gradient(const Point<2> &p, const unsigned int component) const override;
462 
463  virtual void
464  gradient_list(const std::vector<Point<2>> &points,
465  std::vector<Tensor<1, 2>> & gradients,
466  const unsigned int component) const override;
467 
468  virtual void
470  const std::vector<Point<2>> &,
471  std::vector<std::vector<Tensor<1, 2>>> &) const override;
472 
473  virtual double
474  laplacian(const Point<2> &p, const unsigned int component) const override;
475 
476  virtual void
477  laplacian_list(const std::vector<Point<2>> &points,
478  std::vector<double> & values,
479  const unsigned int component) const override;
480  };
481 
482 
483 
489  template <int dim>
490  class SlitSingularityFunction : public Function<dim>
491  {
492  public:
493  virtual double
494  value(const Point<dim> &p, const unsigned int component = 0) const override;
495 
496  virtual void
497  value_list(const std::vector<Point<dim>> &points,
498  std::vector<double> & values,
499  const unsigned int component = 0) const override;
500 
501  virtual void
502  vector_value_list(const std::vector<Point<dim>> &points,
503  std::vector<Vector<double>> & values) const override;
504 
505  virtual Tensor<1, dim>
506  gradient(const Point<dim> & p,
507  const unsigned int component = 0) const override;
508 
509  virtual void
510  gradient_list(const std::vector<Point<dim>> &points,
511  std::vector<Tensor<1, dim>> & gradients,
512  const unsigned int component = 0) const override;
513 
514  virtual void
516  const std::vector<Point<dim>> &,
517  std::vector<std::vector<Tensor<1, dim>>> &) const override;
518 
519  virtual double
520  laplacian(const Point<dim> & p,
521  const unsigned int component = 0) const override;
522 
523  virtual void
524  laplacian_list(const std::vector<Point<dim>> &points,
525  std::vector<double> & values,
526  const unsigned int component = 0) const override;
527  };
528 
529 
536  {
537  public:
538  virtual double
539  value(const Point<2> &p, const unsigned int component = 0) const override;
540 
541  virtual void
542  value_list(const std::vector<Point<2>> &points,
543  std::vector<double> & values,
544  const unsigned int component = 0) const override;
545 
546  virtual void
547  vector_value_list(const std::vector<Point<2>> &points,
548  std::vector<Vector<double>> &values) const override;
549 
550  virtual Tensor<1, 2>
551  gradient(const Point<2> & p,
552  const unsigned int component = 0) const override;
553 
554  virtual void
555  gradient_list(const std::vector<Point<2>> &points,
556  std::vector<Tensor<1, 2>> & gradients,
557  const unsigned int component = 0) const override;
558 
559  virtual void
561  const std::vector<Point<2>> &,
562  std::vector<std::vector<Tensor<1, 2>>> &) const override;
563 
564  virtual double
565  laplacian(const Point<2> & p,
566  const unsigned int component = 0) const override;
567 
568  virtual void
569  laplacian_list(const std::vector<Point<2>> &points,
570  std::vector<double> & values,
571  const unsigned int component = 0) const override;
572  };
573 
574 
575 
590  template <int dim>
591  class JumpFunction : public Function<dim>
592  {
593  public:
598  JumpFunction(const Point<dim> &direction, const double steepness);
599 
603  virtual double
604  value(const Point<dim> &p, const unsigned int component = 0) const override;
605 
609  virtual void
610  value_list(const std::vector<Point<dim>> &points,
611  std::vector<double> & values,
612  const unsigned int component = 0) const override;
613 
617  virtual Tensor<1, dim>
618  gradient(const Point<dim> & p,
619  const unsigned int component = 0) const override;
620 
624  virtual void
625  gradient_list(const std::vector<Point<dim>> &points,
626  std::vector<Tensor<1, dim>> & gradients,
627  const unsigned int component = 0) const override;
628 
632  virtual double
633  laplacian(const Point<dim> & p,
634  const unsigned int component = 0) const override;
635 
639  virtual void
640  laplacian_list(const std::vector<Point<dim>> &points,
641  std::vector<double> & values,
642  const unsigned int component = 0) const override;
643 
650  std::size_t
651  memory_consumption() const;
652 
653  protected:
658 
662  const double steepness;
663 
667  double angle;
668 
672  double sine;
673 
677  double cosine;
678  };
679 
680 
681 
693  template <int dim>
694  class FourierCosineFunction : public Function<dim>
695  {
696  public:
701  FourierCosineFunction(const Tensor<1, dim> &fourier_coefficients);
702 
709  virtual double
710  value(const Point<dim> &p, const unsigned int component = 0) const override;
711 
716  virtual Tensor<1, dim>
717  gradient(const Point<dim> & p,
718  const unsigned int component = 0) const override;
719 
723  virtual double
724  laplacian(const Point<dim> & p,
725  const unsigned int component = 0) const override;
726 
727  private:
732  };
733 
734 
735 
747  template <int dim>
748  class FourierSineFunction : public Function<dim>
749  {
750  public:
755  FourierSineFunction(const Tensor<1, dim> &fourier_coefficients);
756 
763  virtual double
764  value(const Point<dim> &p, const unsigned int component = 0) const override;
765 
770  virtual Tensor<1, dim>
771  gradient(const Point<dim> & p,
772  const unsigned int component = 0) const override;
773 
777  virtual double
778  laplacian(const Point<dim> & p,
779  const unsigned int component = 0) const override;
780 
781  private:
786  };
787 
788 
797  template <int dim>
798  class FourierSineSum : public Function<dim>
799  {
800  public:
805  FourierSineSum(const std::vector<Point<dim>> &fourier_coefficients,
806  const std::vector<double> & weights);
807 
814  virtual double
815  value(const Point<dim> &p, const unsigned int component = 0) const override;
816 
821  virtual Tensor<1, dim>
822  gradient(const Point<dim> & p,
823  const unsigned int component = 0) const override;
824 
828  virtual double
829  laplacian(const Point<dim> & p,
830  const unsigned int component = 0) const override;
831 
832  private:
836  const std::vector<Point<dim>> fourier_coefficients;
837  const std::vector<double> weights;
838  };
839 
840 
841 
851  template <int dim>
852  class FourierCosineSum : public Function<dim>
853  {
854  public:
859  FourierCosineSum(const std::vector<Point<dim>> &fourier_coefficients,
860  const std::vector<double> & weights);
861 
868  virtual double
869  value(const Point<dim> &p, const unsigned int component = 0) const override;
870 
875  virtual Tensor<1, dim>
876  gradient(const Point<dim> & p,
877  const unsigned int component = 0) const override;
878 
882  virtual double
883  laplacian(const Point<dim> & p,
884  const unsigned int component = 0) const override;
885 
886  private:
890  const std::vector<Point<dim>> fourier_coefficients;
891  const std::vector<double> weights;
892  };
893 
894 
906  template <int dim>
907  class CutOffFunctionBase : public Function<dim>
908  {
909  public:
914  static const unsigned int no_component = numbers::invalid_unsigned_int;
915 
932  const double radius = 1.,
933  const Point<dim> center = Point<dim>(),
934  const unsigned int n_components = 1,
935  const unsigned int select = CutOffFunctionBase<dim>::no_component,
936  const bool integrate_to_one = false,
937  const double unitary_integral_value = 1.0);
938 
942  virtual ~CutOffFunctionBase() = default;
943 
947  virtual void
948  set_center(const Point<dim> &p);
949 
953  virtual void
954  set_radius(const double r);
955 
959  const Point<dim> &
960  get_center() const;
961 
965  double
966  get_radius() const;
967 
971  bool
972  integrates_to_one() const;
973 
974  protected:
979 
983  double radius;
984 
989  const unsigned int selected;
990 
995 
1001 
1005  double rescaling;
1006  };
1007 
1008 
1018  template <int dim>
1020  {
1021  public:
1032  double radius = 1.0,
1033  const Point<dim> & center = Point<dim>(),
1034  const unsigned int n_components = 1,
1035  const unsigned int select = CutOffFunctionBase<dim>::no_component,
1036  const bool integrate_to_one = false);
1037 
1042  template <template <int> class CutOffFunctionBaseType>
1043  void
1044  set_base();
1045 
1049  virtual void
1050  set_center(const Point<dim> &center) override;
1051 
1055  virtual void
1056  set_radius(const double radius) override;
1057 
1061  virtual double
1062  value(const Point<dim> &p, const unsigned int component = 0) const override;
1063 
1067  virtual Tensor<1, dim>
1068  gradient(const Point<dim> & p,
1069  const unsigned int component = 0) const override;
1070 
1071  private:
1072  std::array<std::unique_ptr<CutOffFunctionBase<1>>, dim> base;
1073 
1075  };
1076 
1077 
1078 
1087  template <int dim>
1089  {
1090  public:
1098  const double radius = 1.,
1099  const Point<dim> = Point<dim>(),
1100  const unsigned int n_components = 1,
1101  const unsigned int select = CutOffFunctionBase<dim>::no_component,
1102  const bool integrate_to_one = false);
1103 
1107  virtual double
1108  value(const Point<dim> &p, const unsigned int component = 0) const override;
1109 
1113  virtual void
1114  value_list(const std::vector<Point<dim>> &points,
1115  std::vector<double> & values,
1116  const unsigned int component = 0) const override;
1117 
1121  virtual void
1122  vector_value_list(const std::vector<Point<dim>> &points,
1123  std::vector<Vector<double>> & values) const override;
1124  };
1125 
1126 
1135  template <int dim>
1137  {
1138  public:
1146  const double radius = 1.,
1147  const Point<dim> = Point<dim>(),
1148  const unsigned int n_components = 1,
1149  const unsigned int select = CutOffFunctionBase<dim>::no_component,
1150  const bool integrate_to_one = false);
1151 
1155  virtual double
1156  value(const Point<dim> &p, const unsigned int component = 0) const override;
1157 
1161  virtual void
1162  value_list(const std::vector<Point<dim>> &points,
1163  std::vector<double> & values,
1164  const unsigned int component = 0) const override;
1165 
1169  virtual void
1170  vector_value_list(const std::vector<Point<dim>> &points,
1171  std::vector<Vector<double>> & values) const override;
1172  };
1173 
1174 
1187  template <int dim>
1189  {
1190  public:
1195  const double radius = 1.,
1196  const Point<dim> = Point<dim>(),
1197  const unsigned int n_components = 1,
1198  const unsigned int select = CutOffFunctionBase<dim>::no_component,
1199  bool integrate_to_one = false);
1200 
1204  virtual double
1205  value(const Point<dim> &p, const unsigned int component = 0) const override;
1206 
1210  virtual void
1211  value_list(const std::vector<Point<dim>> &points,
1212  std::vector<double> & values,
1213  const unsigned int component = 0) const override;
1214 
1218  virtual void
1219  vector_value_list(const std::vector<Point<dim>> &points,
1220  std::vector<Vector<double>> & values) const override;
1221 
1225  virtual Tensor<1, dim>
1226  gradient(const Point<dim> & p,
1227  const unsigned int component = 0) const override;
1228  };
1229 
1230 
1240  template <int dim>
1242  {
1243  public:
1251  const double radius = 1.,
1252  const Point<dim> = Point<dim>(),
1253  const unsigned int n_components = 1,
1254  const unsigned int select = CutOffFunctionBase<dim>::no_component,
1255  bool integrate_to_one = false);
1256 
1260  virtual double
1261  value(const Point<dim> &p, const unsigned int component = 0) const override;
1262 
1266  virtual void
1267  value_list(const std::vector<Point<dim>> &points,
1268  std::vector<double> & values,
1269  const unsigned int component = 0) const override;
1270 
1274  virtual void
1275  vector_value_list(const std::vector<Point<dim>> &points,
1276  std::vector<Vector<double>> & values) const override;
1277 
1281  virtual Tensor<1, dim>
1282  gradient(const Point<dim> & p,
1283  const unsigned int component = 0) const override;
1284  };
1285 
1286 
1287 
1301  template <int dim>
1302  class Monomial : public Function<dim>
1303  {
1304  public:
1311  Monomial(const Tensor<1, dim> &exponents,
1312  const unsigned int n_components = 1);
1313 
1317  virtual double
1318  value(const Point<dim> &p, const unsigned int component = 0) const override;
1319 
1326  virtual void
1327  vector_value(const Point<dim> &p, Vector<double> &values) const override;
1328 
1332  virtual void
1333  value_list(const std::vector<Point<dim>> &points,
1334  std::vector<double> & values,
1335  const unsigned int component = 0) const override;
1336 
1340  virtual Tensor<1, dim>
1341  gradient(const Point<dim> & p,
1342  const unsigned int component = 0) const override;
1343 
1344  private:
1349  };
1350 
1351 
1352 
1385  template <int dim>
1387  {
1388  public:
1407  const std::array<std::vector<double>, dim> &coordinate_values,
1408  const Table<dim, double> & data_values);
1409 
1420  virtual double
1421  value(const Point<dim> &p, const unsigned int component = 0) const override;
1422 
1434  virtual Tensor<1, dim>
1435  gradient(const Point<dim> & p,
1436  const unsigned int component = 0) const override;
1437 
1438  protected:
1443  table_index_of_point(const Point<dim> &p) const;
1444 
1448  const std::array<std::vector<double>, dim> coordinate_values;
1449 
1454  };
1455 
1456 
1492  template <int dim>
1494  {
1495  public:
1511  const std::array<std::pair<double, double>, dim> &interval_endpoints,
1512  const std::array<unsigned int, dim> & n_subintervals,
1513  const Table<dim, double> & data_values);
1514 
1525  virtual double
1526  value(const Point<dim> &p, const unsigned int component = 0) const override;
1527 
1539  virtual Tensor<1, dim>
1540  gradient(const Point<dim> & p,
1541  const unsigned int component = 0) const override;
1542 
1543  private:
1547  const std::array<std::pair<double, double>, dim> interval_endpoints;
1548 
1552  const std::array<unsigned int, dim> n_subintervals;
1553 
1558  };
1559 
1560 
1573  template <int dim>
1574  class Polynomial : public Function<dim>
1575  {
1576  public:
1586  Polynomial(const Table<2, double> & exponents,
1587  const std::vector<double> &coefficients);
1588 
1592  virtual double
1593  value(const Point<dim> &p, const unsigned int component = 0) const override;
1594 
1595 
1599  virtual void
1600  value_list(const std::vector<Point<dim>> &points,
1601  std::vector<double> & values,
1602  const unsigned int component = 0) const override;
1603 
1607  virtual Tensor<1, dim>
1608  gradient(const Point<dim> & p,
1609  const unsigned int component = 0) const override;
1610 
1611  private:
1616 
1620  const std::vector<double> coefficients;
1621  };
1622 
1623 #ifndef DOXYGEN
1624 
1625 
1626 
1627  // Template definitions
1628  template <int dim>
1629  template <template <int> class CutOffFunctionBaseType>
1630  void
1632  {
1633  initialized = true;
1634  static_assert(
1635  std::is_base_of<CutOffFunctionBase<1>, CutOffFunctionBaseType<1>>::value,
1636  "You can only construct a CutOffFunctionTensorProduct function from "
1637  "a class derived from CutOffFunctionBase.");
1638 
1639  for (unsigned int i = 0; i < dim; ++i)
1640  base[i].reset(new CutOffFunctionBaseType<1>(this->radius,
1641  Point<1>(this->center[i]),
1642  this->n_components,
1643  this->selected,
1644  this->integrate_to_one));
1645  }
1646 
1647 
1648 
1649 #endif
1650 
1651 } // namespace Functions
1653 
1654 #endif
const Tensor< 1, dim > exponents
const std::vector< double > weights
Definition: function_lib.h:837
static const unsigned int invalid_unsigned_int
Definition: types.h:196
const unsigned int n_components
Definition: function.h:164
virtual void vector_gradient_list(const std::vector< Point< dim >> &points, std::vector< std::vector< Tensor< 1, dim, double >>> &gradients) const
virtual void value_list(const std::vector< Point< dim >> &points, std::vector< double > &values, const unsigned int component=0) const override
Definition: function_lib.cc:52
virtual void vector_gradient(const Point< dim > &p, std::vector< Tensor< 1, dim >> &gradient) const override
std::array< std::unique_ptr< CutOffFunctionBase< 1 > >, dim > base
virtual double value(const Point< dim > &p, const unsigned int component=0) const override
Definition: function_lib.cc:34
virtual void gradient_list(const std::vector< Point< dim >> &points, std::vector< Tensor< 1, dim >> &gradients, const unsigned int component=0) const override
virtual Tensor< 1, dim > gradient(const Point< dim > &p, const unsigned int component=0) const override
Definition: function_lib.cc:92
const Tensor< 1, dim > fourier_coefficients
Definition: function_lib.h:731
virtual void vector_value_list(const std::vector< Point< dim >> &points, std::vector< Vector< double >> &values) const
virtual void vector_value(const Point< dim > &p, Vector< double > &values) const override
Definition: function_lib.cc:42
virtual SymmetricTensor< 2, dim, double > hessian(const Point< dim > &p, const unsigned int component=0) const
const Table< dim, double > data_values
const std::vector< double > weights
Definition: function_lib.h:891
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
const std::vector< Point< dim > > fourier_coefficients
Definition: function_lib.h:890
virtual void hessian_list(const std::vector< Point< dim >> &points, std::vector< SymmetricTensor< 2, dim, double >> &values, const unsigned int component=0) const
virtual double laplacian(const Point< dim > &p, const unsigned int component=0) const override
Definition: function_lib.cc:69
const Table< 2, double > exponents
const std::vector< double > coefficients
const unsigned int selected
Definition: function_lib.h:989
const Point< dim > direction
Definition: function_lib.h:657
const std::vector< Point< dim > > fourier_coefficients
Definition: function_lib.h:836
const std::array< std::vector< double >, dim > coordinate_values
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
const std::array< unsigned int, dim > n_subintervals
Point< 3 > center
const Tensor< 1, dim > fourier_coefficients
Definition: function_lib.h:785
std::size_t memory_consumption() const
const std::array< std::pair< double, double >, dim > interval_endpoints
virtual void laplacian_list(const std::vector< Point< dim >> &points, std::vector< double > &values, const unsigned int component=0) const override
Definition: function_lib.cc:77