v0.14.0
BiLinearFormsIntegratorsImpl.hpp
Go to the documentation of this file.
1 /** \file BiLinearFormsIntegratorsImpl.hpp
2  * \brief Bilinear forms integrators (implementation)
3  * \ingroup mofem_form
4 
5  \todo Some operators could be optimised. To do that, we need to write tests
6  and use Valgrind to profile code, shaking cache misses. For example, some
7  operators should have iteration first over columns, then rows. ome operators.
8  Since those operators are used in many problems, an implementation must be
9  efficient.
10 
11 */
12 
13 #ifndef __BILINEAR_FORMS_INTEGRATORS_IMPL_HPP__
14 #define __BILINEAR_FORMS_INTEGRATORS_IMPL_HPP__
15 
16 namespace MoFEM {
17 
18 using CacheMatsTypeType =
19  std::map<std::pair<int, int>, boost::shared_ptr<MatrixDouble>>;
20 
21 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
22  typename OpBase>
23 struct OpGradGradImpl {};
24 
25 template <int SPACE_DIM, typename OpBase>
26 struct OpGradGradImpl<1, 1, SPACE_DIM, GAUSS, OpBase> : public OpBase {
27  FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
28  OpGradGradImpl(const std::string row_field_name,
29  const std::string col_field_name,
31  boost::shared_ptr<Range> ents_ptr = nullptr)
32  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
33  betaCoeff(beta) {
34  if (row_field_name == col_field_name)
35  this->sYmm = true;
36  }
37 
38 protected:
40  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
41  EntitiesFieldData::EntData &col_data);
42 };
43 
44 template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
46  FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
47  OpGradGradImpl(const std::string row_field_name,
48  const std::string col_field_name,
50  boost::shared_ptr<Range> ents_ptr = nullptr)
51  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
52  betaCoeff(beta) {
53  if (row_field_name == col_field_name)
54  this->sYmm = true;
55  }
56 
57 protected:
59  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
60  EntitiesFieldData::EntData &col_data);
61 };
62 
63 template <int BASE_DIM, int FIELD_DIM, IntegrationType I, typename OpBase>
64 struct OpMassImpl {};
65 
66 template <typename OpBase>
67 struct OpMassImpl<1, 1, GAUSS, OpBase> : public OpBase {
68 
70  const std::string row_field_name, const std::string col_field_name,
71  ScalarFun beta = [](double, double, double) constexpr { return 1; },
72  boost::shared_ptr<Range> ents_ptr = nullptr,
73  boost::shared_ptr<MatrixDouble> cache_mat = nullptr)
74  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
75  betaCoeff(beta) {
76  if (row_field_name == col_field_name)
77  this->sYmm = true;
78  }
79 
80 protected:
81  ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
82  MoFEMErrorCode integrateImpl(EntitiesFieldData::EntData &row_data,
84  double vol);
86  EntitiesFieldData::EntData &col_data) {
87  return integrateImpl(row_data, col_data, OpBase::getMeasure());
88  }
89 };
90 
91 template <int FIELD_DIM, typename OpBase>
93  : public OpMassImpl<1, 1, GAUSS, OpBase> {
95 
96 protected:
98  MoFEMErrorCode integrateImpl(EntitiesFieldData::EntData &row_data,
100  double vol);
102  EntitiesFieldData::EntData &col_data) {
103  return integrateImpl(row_data, col_data, OpBase::getMeasure());
104  }
105 };
106 
107 template <int FIELD_DIM, typename OpBase>
108 struct OpMassImpl<3, FIELD_DIM, GAUSS, OpBase> : public OpBase {
109 
110  OpMassImpl(const std::string row_field_name, const std::string col_field_name,
111  ScalarFun beta = scalar_fun_one,
112  boost::shared_ptr<Range> ents_ptr = nullptr)
113  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
114  betaCoeff(beta) {
115  if (row_field_name == col_field_name)
116  this->sYmm = true;
117  }
118 
119 protected:
121  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
122  EntitiesFieldData::EntData &col_data);
123 };
124 
125 template <typename OpBase>
126 struct OpMassImpl<3, 4, GAUSS, OpBase> : public OpBase {
127  OpMassImpl(const std::string row_field_name, const std::string col_field_name,
128  ScalarFun beta = scalar_fun_one,
129  boost::shared_ptr<Range> ents_ptr = nullptr)
130  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
131  betaCoeff(beta) {
132  if (row_field_name == col_field_name)
133  this->sYmm = true;
134  }
135 
136 protected:
138  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
139  EntitiesFieldData::EntData &col_data);
140 };
141 
142 template <typename OpBase>
143 struct OpMassImpl<3, 9, GAUSS, OpBase> : public OpBase {
144  OpMassImpl(const std::string row_field_name, const std::string col_field_name,
145  ScalarFun beta = scalar_fun_one,
146  boost::shared_ptr<Range> ents_ptr = nullptr)
147  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
148  betaCoeff(beta) {
149  if (row_field_name == col_field_name)
150  this->sYmm = true;
151  }
152 
153 protected:
155  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
156  EntitiesFieldData::EntData &col_data);
157 };
158 
159 template <int BASE_DIM, int FIELD_DIM, IntegrationType I, typename OpBase>
161 
162 template <int FIELD_DIM, IntegrationType I, typename OpBase>
164  : public OpMassImpl<1, FIELD_DIM, I, OpBase> {
165  OpMassCacheImpl(CacheMatsTypeType cache, const std::string row_field_name,
166  const std::string col_field_name, const double beta,
167  boost::shared_ptr<Range> ents_ptr = nullptr)
168  : OpMassImpl<1, FIELD_DIM, I, OpBase>(
169  row_field_name, col_field_name,
170  [](double, double, double) constexpr { return 1; }, ents_ptr),
171  cacheLocMats(cache), scalarBeta(beta) {}
172 
173  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
174  EntitiesFieldData::EntData &col_data);
175 
176 protected:
177  double scalarBeta;
179 };
180 
181 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
182  typename OpBase>
184 
185 template <int SPACE_DIM, int S, typename OpBase>
187  : public OpBase {
189  OpGradSymTensorGradImpl(const std::string row_field_name,
190  const std::string col_field_name,
191  boost::shared_ptr<MatrixDouble> mat_D,
192  boost::shared_ptr<Range> ents_ptr = nullptr,
193  ScalarFun beta = scalar_fun_one)
194  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
195  matD(mat_D), betaCoeff(beta) {
196  if (row_field_name == col_field_name)
197  this->sYmm = true;
198  }
199 
200 protected:
201  boost::shared_ptr<MatrixDouble> matD;
203  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
204  EntitiesFieldData::EntData &col_data);
205 };
206 
207 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
208  typename OpBase>
210 
211 template <int SPACE_DIM, int S, typename OpBase>
213  : public OpBase {
215  OpGradTensorGradImpl(const std::string row_field_name,
216  const std::string col_field_name,
217  boost::shared_ptr<MatrixDouble> mat_D,
218  ScalarFun beta = scalar_fun_one)
219  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL), matD(mat_D),
220  betaCoeff(beta) {}
221 
222 protected:
223  boost::shared_ptr<MatrixDouble> matD;
225  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
226  EntitiesFieldData::EntData &col_data);
227 };
228 
229 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
230  typename OpBase>
232 
233 template <int SPACE_DIM, int S, typename OpBase>
235  : public OpBase {
237  OpGradGradSymTensorGradGradImpl(const std::string row_field_name,
238  const std::string col_field_name,
239  boost::shared_ptr<MatrixDouble> mat_D,
240  boost::shared_ptr<Range> ents_ptr = nullptr)
241  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
242  matD(mat_D) {
243  if (row_field_name == col_field_name)
244  this->sYmm = true;
245  }
246 
247 protected:
248  boost::shared_ptr<MatrixDouble> matD;
249  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
250  EntitiesFieldData::EntData &col_data);
251 };
252 
253 template <int SPACE_DIM, IntegrationType I, typename OpBase>
255 
256 template <int SPACE_DIM, typename OpBase>
258  OpMixDivTimesScalarImpl(const std::string row_field_name,
259  const std::string col_field_name,
260  ConstantFun alpha_fun,
261  const bool assemble_transpose = false,
262  const bool only_transpose = false)
263  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
264  alphaConstant(alpha_fun) {
265  this->assembleTranspose = assemble_transpose;
266  this->onlyTranspose = only_transpose;
267  }
268 
269 protected:
272  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
273  EntitiesFieldData::EntData &col_data);
274 };
275 
276 template <int SPACE_DIM, IntegrationType I, typename OpBase,
277  CoordinateTypes CoordSys>
279 
280 template <int SPACE_DIM, typename OpBase, CoordinateTypes CoordSys>
282  : public OpBase {
283 
284  OpMixDivTimesVecImpl(const std::string row_field_name,
285  const std::string col_field_name, ConstantFun alpha_fun,
286  const bool assemble_transpose,
287  const bool only_transpose = false)
288  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
289  alphaConstant(alpha_fun) {
290  this->assembleTranspose = assemble_transpose;
291  this->onlyTranspose = only_transpose;
292  }
293 
294  OpMixDivTimesVecImpl(const std::string row_field_name,
295  const std::string col_field_name, ConstantFun alpha_fun,
296  ScalarFun beta_fun, const bool assemble_transpose,
297  const bool only_transpose = false)
298  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
299  alphaConstant(alpha_fun), betaCoeff(beta_fun) {
300  this->assembleTranspose = assemble_transpose;
301  this->onlyTranspose = only_transpose;
302  }
303 
304 protected:
306 
307  ConstantFun alphaConstant = []() constexpr { return 1; };
308  ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
309 
310  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
311  EntitiesFieldData::EntData &col_data);
312 };
313 
314 template <int SPACE_DIM, IntegrationType I, typename OpBase,
315  CoordinateTypes COORDINATE_SYSTEM>
317 
318 template <int SPACE_DIM, typename OpBase, CoordinateTypes COORDINATE_SYSTEM>
319 struct OpMixScalarTimesDivImpl<SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM>
320  : public OpBase {
322  const std::string row_field_name, const std::string col_field_name,
323  ScalarFun alpha_fun = [](double, double, double) constexpr { return 1; },
324  const bool assemble_transpose = false, const bool only_transpose = false)
325  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
326  alphaConstant(alpha_fun) {
327  this->assembleTranspose = assemble_transpose;
328  this->onlyTranspose = only_transpose;
329  this->sYmm = false;
330  }
331 
332 protected:
335  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
336  EntitiesFieldData::EntData &col_data);
337 };
338 
339 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
340  typename OpBase>
342 
343 template <int SPACE_DIM, typename OpBase>
345  : public OpBase {
346  OpMixVectorTimesGradImpl(const std::string row_field_name,
347  const std::string col_field_name,
348  ConstantFun alpha_fun,
349  const bool assemble_transpose = false,
350  const bool only_transpose = false)
351  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
352  alphaConstant(alpha_fun) {
353  this->assembleTranspose = assemble_transpose;
354  this->onlyTranspose = only_transpose;
355  }
356 
357 protected:
360  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
361  EntitiesFieldData::EntData &col_data);
362 };
363 
364 template <int SPACE_DIM, typename OpBase>
366  : public OpBase {
367  OpMixVectorTimesGradImpl(const std::string row_field_name,
368  const std::string col_field_name,
369  ConstantFun alpha_fun,
370  const bool assemble_transpose = false,
371  const bool only_transpose = false)
372  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
373  alphaConstant(alpha_fun) {
374  this->assembleTranspose = assemble_transpose;
375  this->onlyTranspose = only_transpose;
376  }
377 
378 protected:
381  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
382  EntitiesFieldData::EntData &col_data);
383 };
384 
385 template <int SPACE_DIM, IntegrationType I, typename OpBase>
387 
388 template <int SPACE_DIM, typename OpBase>
390 
391  OpMixTensorTimesGradImpl(const std::string row_field_name,
392  const std::string col_field_name,
393  ConstantFun alpha_fun, const bool assemble_transpose,
394  const bool only_transpose = false)
395  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
396  alphaConstant(alpha_fun) {
397  this->assembleTranspose = assemble_transpose;
398  this->onlyTranspose = only_transpose;
399  }
400 
401  OpMixTensorTimesGradImpl(const std::string row_field_name,
402  const std::string col_field_name,
403  ConstantFun alpha_fun, ScalarFun beta_coeff,
404  const bool assemble_transpose,
405  const bool only_transpose = false)
406  : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
407  alphaConstant(alpha_fun), betaCoeff(beta_coeff) {
408  this->assembleTranspose = assemble_transpose;
409  this->onlyTranspose = only_transpose;
410  }
411 
412 protected:
415  ConstantFun alphaConstant = []() constexpr { return 1; };
416  ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
417  MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
418  EntitiesFieldData::EntData &col_data);
419 };
420 
421 template <int SPACE_DIM, typename OpBase>
423  EntitiesFieldData::EntData &row_data,
424  EntitiesFieldData::EntData &col_data) {
426 
427  // get element volume
428  const double vol = OpBase::getMeasure();
429  // get integration weights
430  auto t_w = OpBase::getFTensor0IntegrationWeight();
431  // get base function gradient on rows
432  auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
433  // get coordinate at integration points
434  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
435 
436  // loop over integration points
437  for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
438  const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
439  // take into account Jacobian
440  const double alpha = t_w * beta;
441  // loop over ros base functions
442  int rr = 0;
443  for (; rr != OpBase::nbRows; rr++) {
444  // get column base functions gradient at gauss point gg
445  auto t_col_grad = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
446  // loop over columns
447  for (int cc = 0; cc != OpBase::nbCols; cc++) {
448  // calculate element of local matrix
449  OpBase::locMat(rr, cc) += alpha * (t_row_grad(i) * t_col_grad(i));
450  ++t_col_grad; // move to another gradient of base function
451  // on column
452  }
453  ++t_row_grad; // move to another element of gradient of base
454  // function on row
455  }
456  for (; rr < OpBase::nbRowBaseFunctions; ++rr)
457  ++t_row_grad;
458 
459  ++t_coords;
460  ++t_w; // move to another integration weight
461  }
463 }
464 
465 template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
468  EntitiesFieldData::EntData &row_data,
469  EntitiesFieldData::EntData &col_data) {
471 
472  // get element volume
473  const double vol = OpBase::getMeasure();
474  // get integration weights
475  auto t_w = OpBase::getFTensor0IntegrationWeight();
476  // get base function gradient on rows
477  auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
478  // get coordinate at integration points
479  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
480 
483 
484  auto get_t_vec = [&](const int rr) {
485  std::array<double *, FIELD_DIM> ptrs;
486  for (auto i = 0; i != FIELD_DIM; ++i)
487  ptrs[i] = &OpBase::locMat(rr + i, i);
489  ptrs);
490  };
491 
492  // loop over integration points
493  for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
494  const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
495  // take into account Jacobian
496  const double alpha = t_w * beta;
497  // loop over ros base functions
498  int rr = 0;
499  for (; rr != OpBase::nbRows / FIELD_DIM; rr++) {
500  // get diag vec
501  auto t_vec = get_t_vec(rr * FIELD_DIM);
502  // get column base functions gradient at gauss point gg
503  auto t_col_grad = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
504  // loop over columns
505  for (int cc = 0; cc != OpBase::nbCols / FIELD_DIM; cc++) {
506  // calculate element of local matrix
507  t_vec(j) += alpha * (t_row_grad(i) * t_col_grad(i));
508  ++t_col_grad; // move to another gradient of base function
509  // on column
510  ++t_vec;
511  }
512  ++t_row_grad; // move to another element of gradient of base
513  // function on row
514  }
515  for (; rr < OpBase::nbRowBaseFunctions; ++rr)
516  ++t_row_grad;
517 
518  ++t_coords;
519  ++t_w; // move to another integration weight
520  }
522 }
523 
524 template <typename OpBase>
526  EntitiesFieldData::EntData &row_data,
527  EntitiesFieldData::EntData &col_data, double vol) {
529 
530 #ifndef NDEBUG
531  auto log_error = [&]() {
532  MOFEM_LOG("SELF", Sev::error) << "Row side " << OpBase::rowSide << " "
533  << CN::EntityTypeName(OpBase::rowType);
534  MOFEM_LOG("SELF", Sev::error) << "Col side " << OpBase::colSide << " "
535  << CN::EntityTypeName(OpBase::colType);
536  };
537 
538  if (row_data.getN().size2() < OpBase::nbRows) {
539  log_error();
540  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
541  "Wrong number of base functions on rows %d < %d",
542  row_data.getN().size2(), OpBase::nbRows);
543  }
544  if (col_data.getN().size2() < OpBase::nbCols) {
545  log_error();
546  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
547  "Wrong number of base functions on cols %d < %d",
548  col_data.getN().size2(), OpBase::nbCols);
549  }
550  if (row_data.getN().size1() != OpBase::nbIntegrationPts) {
551  log_error();
552  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
553  "Wrong number of integration points on rows %d != %d",
554  row_data.getN().size1(), OpBase::nbIntegrationPts);
555  }
556  if (col_data.getN().size1() != OpBase::nbIntegrationPts) {
557  log_error();
558  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
559  "Wrong number of integration points on cols %d != %d",
560  col_data.getN().size1(), OpBase::nbIntegrationPts);
561  }
562 #endif
563 
564  // get integration weights
565  auto t_w = OpBase::getFTensor0IntegrationWeight();
566  // get base function gradient on rows
567  auto t_row_base = row_data.getFTensor0N();
568  // get coordinate at integration points
569  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
570  // loop over integration points
571  for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
572  const double beta = betaCoeff(t_coords(0), t_coords(1), t_coords(2));
573  // take into account Jacobian
574  const double alpha = t_w * beta;
575  // loop over rows base functions
576  auto a_mat_ptr = &*OpBase::locMat.data().begin();
577  int rr = 0;
578  for (; rr != OpBase::nbRows; rr++) {
579  // get column base functions gradient at gauss point gg
580  auto t_col_base = col_data.getFTensor0N(gg, 0);
581  // loop over columns
582  for (int cc = 0; cc != OpBase::nbCols; cc++) {
583  // calculate element of local matrix
584  *a_mat_ptr += alpha * (t_row_base * t_col_base);
585  ++t_col_base;
586  ++a_mat_ptr;
587  }
588  ++t_row_base;
589  }
590  for (; rr < OpBase::nbRowBaseFunctions; ++rr)
591  ++t_row_base;
592  ++t_coords;
593  ++t_w; // move to another integration weight
594  }
595 
596  OpBase::locMat *= vol;
597 
599 };
600 
601 template <int FIELD_DIM, typename OpBase>
604  double vol) {
606 
607  auto integrate = [&](auto &mat) {
609  // get integration weights
610  auto t_w = OpBase::getFTensor0IntegrationWeight();
611  // get base function gradient on rows
612  auto t_row_base = row_data.getFTensor0N();
613  // get coordinate at integration points
614  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
615 
617  auto get_t_vec = [&](const int rr) {
618  std::array<double *, FIELD_DIM> ptrs;
619  for (auto i = 0; i != FIELD_DIM; ++i)
620  ptrs[i] = &mat(rr + i, i);
622  ptrs);
623  };
624 
625  // loop over integration points
626  for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
627  const double beta = betaCoeff(t_coords(0), t_coords(1), t_coords(2));
628  // take into account Jacobian
629  const double alpha = t_w * beta;
630  // loop over rows base functions
631  int rr = 0;
632  for (; rr != OpBase::nbRows / FIELD_DIM; rr++) {
633  // get column base functions gradient at gauss point gg
634  auto t_col_base = col_data.getFTensor0N(gg, 0);
635  // get mat vec
636  auto t_vec = get_t_vec(FIELD_DIM * rr);
637  // loop over columns
638  for (int cc = 0; cc != OpBase::nbCols / FIELD_DIM; cc++) {
639  // calculate element of local matrix
640  t_vec(i) += alpha * (t_row_base * t_col_base);
641  ++t_col_base;
642  ++t_vec;
643  }
644  ++t_row_base;
645  }
646  for (; rr < OpBase::nbRowBaseFunctions; ++rr)
647  ++t_row_base;
648  ++t_coords;
649  ++t_w; // move to another integration weight
650  }
652  };
653 
654  CHKERR integrate(OpBase::locMat);
655  OpBase::locMat *= vol;
656 
658 };
659 
660 template <int FIELD_DIM, typename OpBase>
662  EntitiesFieldData::EntData &row_data,
663  EntitiesFieldData::EntData &col_data) {
666  size_t nb_base_functions = row_data.getN().size2() / 3;
667  // // get element volume
668  const double vol = OpBase::getMeasure();
669  // get integration weights
670  auto t_w = OpBase::getFTensor0IntegrationWeight();
671  // get base function gradient on rows
672  auto t_row_base = row_data.getFTensor1N<3>();
673  // get coordinate at integration points
674  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
675  // loop over integration points
676  for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
677  const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
678  // take into account Jacobian
679  const double alpha = t_w * beta;
680  // loop over rows base functions
681  auto a_mat_ptr = &*OpBase::locMat.data().begin();
682  int rr = 0;
683  for (; rr != OpBase::nbRows; rr++) {
684  // get column base functions gradient at gauss point gg
685  auto t_col_base = col_data.getFTensor1N<3>(gg, 0);
686  // loop over columns
687  for (int cc = 0; cc != OpBase::nbCols; cc++) {
688  // calculate element of local matrix
689  (*a_mat_ptr) += alpha * (t_row_base(i) * t_col_base(i));
690  ++t_col_base;
691  ++a_mat_ptr;
692  }
693  ++t_row_base;
694  }
695  for (; rr < nb_base_functions; ++rr)
696  ++t_row_base;
697  ++t_coords;
698  ++t_w; // move to another integration weight
699  }
701 };
702 
703 template <typename OpBase>
705  EntitiesFieldData::EntData &row_data,
706  EntitiesFieldData::EntData &col_data) {
710  auto get_t_vec = [&](const int rr) {
712  &OpBase::locMat(rr + 0, 0), &OpBase::locMat(rr + 1, 1)};
713  };
714  size_t nb_base_functions = row_data.getN().size2() / 3;
715  // // get element volume
716  const double vol = OpBase::getMeasure();
717  // get integration weights
718  auto t_w = OpBase::getFTensor0IntegrationWeight();
719  // get base function gradient on rows
720  auto t_row_base = row_data.getFTensor1N<3>();
721  // get coordinate at integration points
722  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
723  // loop over integration points
724  for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
725  const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
726  // take into account Jacobian
727  const double alpha = t_w * beta;
728  // loop over rows base functions
729  int rr = 0;
730  for (; rr != OpBase::nbRows / 2; rr++) {
731  // get column base functions gradient at gauss point gg
732  auto t_col_base = col_data.getFTensor1N<3>(gg, 0);
733  auto t_vec = get_t_vec(2 * rr);
734  // loop over columns
735  for (int cc = 0; cc != OpBase::nbCols / 2; cc++) {
736  // calculate element of local matrix
737  t_vec(I) += alpha * (t_row_base(k) * t_col_base(k));
738  ++t_col_base;
739  ++t_vec;
740  }
741  ++t_row_base;
742  }
743  for (; rr < nb_base_functions; ++rr)
744  ++t_row_base;
745  ++t_coords;
746  ++t_w; // move to another integration weight
747  }
749 }
750 
751 template <typename OpBase>
753  EntitiesFieldData::EntData &row_data,
754  EntitiesFieldData::EntData &col_data) {
758  auto get_t_vec = [&](const int rr) {
760  &OpBase::locMat(rr + 0, 0), &OpBase::locMat(rr + 1, 1),
761  &OpBase::locMat(rr + 2, 2)};
762  };
763  size_t nb_base_functions = row_data.getN().size2() / 3;
764  // // get element volume
765  const double vol = OpBase::getMeasure();
766  // get integration weights
767  auto t_w = OpBase::getFTensor0IntegrationWeight();
768  // get base function gradient on rows
769  auto t_row_base = row_data.getFTensor1N<3>();
770  // get coordinate at integration points
771  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
772  // loop over integration points
773  for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
774  const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
775  // take into account Jacobian
776  const double alpha = t_w * beta;
777  // loop over rows base functions
778  int rr = 0;
779  for (; rr != OpBase::nbRows / 3; rr++) {
780  // get column base functions gradient at gauss point gg
781  auto t_col_base = col_data.getFTensor1N<3>(gg, 0);
782  auto t_vec = get_t_vec(3 * rr);
783  // loop over columns
784  for (int cc = 0; cc != OpBase::nbCols / 3; cc++) {
785  // calculate element of local matrix
786  t_vec(i) += alpha * (t_row_base(k) * t_col_base(k));
787  ++t_col_base;
788  ++t_vec;
789  }
790  ++t_row_base;
791  }
792  for (; rr < nb_base_functions; ++rr)
793  ++t_row_base;
794  ++t_coords;
795  ++t_w; // move to another integration weight
796  }
798 };
799 
800 template <int FIELD_DIM, IntegrationType I, typename OpBase>
802  EntitiesFieldData::EntData &row_data,
803  EntitiesFieldData::EntData &col_data) {
805  CHKERR this->integrateImpl(row_data, col_data, 1);
806 
807  const auto vol = this->getMeasure();
808  const auto row_type = this->rowType;
809  const auto col_type = this->colType;
810  auto &loc_mat = this->locMat;
811 
812  auto p = std::make_pair(row_type, col_type);
813 
814  if (cacheLocMats[p]) {
815  if (cacheLocMats[p]->size1() != loc_mat.size1() &&
816  cacheLocMats[p]->size2() != loc_mat.size2()) {
817  cacheLocMats[p]->resize(loc_mat.size1(), loc_mat.size2());
818  *(cacheLocMats[p]) = loc_mat;
819  }
820  loc_mat = *(cacheLocMats[p]);
821  loc_mat *= scalarBeta * this->getMeasure();
822  } else {
823  loc_mat *= scalarBeta * this->getMeasure();
824  }
826 }
827 
828 template <int SPACE_DIM, int S, typename OpBase>
831  EntitiesFieldData::EntData &row_data,
832  EntitiesFieldData::EntData &col_data) {
834 
835  const size_t nb_row_dofs = row_data.getIndices().size();
836  const size_t nb_col_dofs = col_data.getIndices().size();
837 
838  if (nb_row_dofs && nb_col_dofs) {
839 
840  const bool diag = (row_data.getFieldEntities()[0]->getLocalUniqueId() ==
841  col_data.getFieldEntities()[0]->getLocalUniqueId());
842 
847 
848  // get element volume
849  double vol = OpBase::getMeasure();
850 
851  // get intergrayion weights
852  auto t_w = OpBase::getFTensor0IntegrationWeight();
853 
854  // get derivatives of base functions on rows
855  auto t_row_diff_base = row_data.getFTensor1DiffN<SPACE_DIM>();
856 
857  // Elastic stiffness tensor (4th rank tensor with minor and major
858  // symmetry)
859  auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, S>(*matD);
860 
861  // get coordinate at integration points
862  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
863 
864  // iterate over integration points
865  for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
866 
867  // calculate scalar weight times element volume
868  double a = t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
869 
870  // iterate over row base functions
871  int rr = 0;
872  for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
873 
874  // get sub matrix for the row
875  auto t_m = OpBase::template getLocMat<SPACE_DIM>(SPACE_DIM * rr);
876 
878  // I mix up the indices here so that it behaves like a
879  // Dg. That way I don't have to have a separate wrapper
880  // class Christof_Expr, which simplifies things.
881  t_rowD(l, j, k) = t_D(i, j, k, l) * (a * t_row_diff_base(i));
882 
883  // get derivatives of base functions for columns
884  auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
885 
886  // iterate column base functions
887  const auto nb_cols = (diag) ? rr : OpBase::nbCols / SPACE_DIM - 1;
888  for (int cc = 0; cc <= nb_cols; ++cc) {
889 
890  // integrate block local stiffens matrix
891  t_m(i, j) += t_rowD(i, j, k) * t_col_diff_base(k);
892 
893  // move to next column base function
894  ++t_col_diff_base;
895 
896  // move to next block of local stiffens matrix
897  ++t_m;
898  }
899 
900  // move to next row base function
901  ++t_row_diff_base;
902  }
903 
904  for (; rr < OpBase::nbRowBaseFunctions; ++rr)
905  ++t_row_diff_base;
906 
907  // move to next integration weight
908  ++t_w;
909  ++t_D;
910  ++t_coords;
911  }
912 
913  // Copy symmetry
914  if (diag) {
915  for (int rr = 0; rr != OpBase::nbRows / SPACE_DIM - 1; ++rr) {
916  auto t_m_rr = getFTensor2FromArray<SPACE_DIM, SPACE_DIM, SPACE_DIM>(
917  this->locMat, SPACE_DIM * rr, SPACE_DIM * (rr + 1));
918  auto t_m_cc = getFTensor2FromArray<SPACE_DIM, SPACE_DIM>(
919  this->locMat, SPACE_DIM * (rr + 1), SPACE_DIM * rr,
921  for (int cc = rr + 1; cc < OpBase::nbCols / SPACE_DIM; ++cc) {
922  t_m_rr(i, j) = t_m_cc(j, i);
923  ++t_m_rr;
924  ++t_m_cc;
925  }
926  }
927  }
928  }
929 
931 }
932 
933 template <int SPACE_DIM, int S, typename OpBase>
936  EntitiesFieldData::EntData &row_data,
937  EntitiesFieldData::EntData &col_data) {
939 
941 
946 
947  auto &row_hessian = row_data.getN(BaseDerivatives::SecondDerivative);
948  auto &col_hessian = col_data.getN(BaseDerivatives::SecondDerivative);
949 
950 #ifndef NDEBUG
951  if (row_hessian.size1() != OpBase::nbIntegrationPts) {
952  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
953  "Wrong number of integration pts (%d != %d)",
954  row_hessian.size1(), OpBase::nbIntegrationPts);
955  }
956  if (row_hessian.size2() !=
958  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
959  "Wrong number of base functions (%d != %d)",
960  row_hessian.size2() / (SPACE_DIM * SPACE_DIM),
962  }
963  if (row_hessian.size2() < OpBase::nbRows * SPACE_DIM * SPACE_DIM) {
964  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
965  "Wrong number of base functions (%d < %d)", row_hessian.size2(),
967  }
968  if (col_hessian.size1() != OpBase::nbIntegrationPts) {
969  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
970  "Wrong number of integration pts (%d != %d)",
971  col_hessian.size1(), OpBase::nbIntegrationPts);
972  }
973  if (col_hessian.size2() < OpBase::nbCols * SPACE_DIM * SPACE_DIM) {
974  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
975  "Wrong number of base functions (%d < %d)", col_hessian.size2(),
977  }
978 #endif
979 
980  // get element volume
981  double vol = OpBase::getMeasure();
982 
983  // get intergrayion weights
984  auto t_w = OpBase::getFTensor0IntegrationWeight();
985 
986  auto t_row_diff2 = getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(
987  &*row_hessian.data().begin());
988 
989  // Elastic stiffness tensor (4th rank tensor with minor and major
990  // symmetry)
991  auto t_D = getFTensor4DdgFromMat<SPACE_DIM, SPACE_DIM, S>(*matD);
992 
993  // iterate over integration points
994  for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
995 
996  // calculate scalar weight times element volume
997  double a = t_w * vol;
998 
999  // get sub matrix for the row
1000  auto m_ptr = &*OpBase::locMat.data().begin();
1001 
1002  // iterate over row base functions
1003  int rr = 0;
1004  for (; rr != OpBase::nbRows; ++rr) {
1005 
1007  t_rowD(k, l) = t_D(i, j, k, l) * (a * t_row_diff2(i, j));
1008 
1009  // get derivatives of base functions for columns
1010  auto t_col_diff2 =
1011  getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(&col_hessian(gg, 0));
1012 
1013  // iterate column base functions
1014  for (int cc = 0; cc != OpBase::nbCols; ++cc) {
1015 
1016  // integrate block local stiffens matrix
1017  *m_ptr += t_rowD(i, j) * t_col_diff2(i, j);
1018 
1019  // move to next column base function
1020  ++t_col_diff2;
1021 
1022  // move to next block of local stiffens matrix
1023  ++m_ptr;
1024  }
1025 
1026  // move to next row base function
1027  ++t_row_diff2;
1028  }
1029 
1030  for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1031  ++t_row_diff2;
1032 
1033  // move to next integration weight
1034  ++t_w;
1035  ++t_D;
1036  }
1037  }
1038 
1040 }
1041 
1042 template <int SPACE_DIM, int S, typename OpBase>
1045  EntitiesFieldData::EntData &row_data,
1046  EntitiesFieldData::EntData &col_data) {
1048 
1049  const size_t nb_row_dofs = row_data.getIndices().size();
1050  const size_t nb_col_dofs = col_data.getIndices().size();
1051 
1052  if (nb_row_dofs && nb_col_dofs) {
1053 
1058 
1059  // get element volume
1060  double vol = OpBase::getMeasure();
1061 
1062  // get intergrayion weights
1063  auto t_w = OpBase::getFTensor0IntegrationWeight();
1064 
1065  // get derivatives of base functions on rows
1066  auto t_row_diff_base = row_data.getFTensor1DiffN<SPACE_DIM>();
1067 
1068  // stiffness tensor (4th rank tensor)
1069  auto t_D =
1070  getFTensor4FromMat<SPACE_DIM, SPACE_DIM, SPACE_DIM, SPACE_DIM, S>(
1071  *matD);
1072 
1073  // get coordinate at integration points
1074  auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1075 
1076  // iterate over integration points
1077  for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1078 
1079  // calculate scalar weight times element volume
1080  double a = t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1081 
1082  // iterate over row base functions
1083  int rr = 0;
1084  for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1085 
1086  // get sub matrix for the row
1087  auto t_m = OpBase::template getLocMat<SPACE_DIM>(SPACE_DIM * rr);
1088 
1089  // calculate row
1091  t_row(i, k, l) = t_D(i, j, k, l) * (a * t_row_diff_base(j));
1092 
1093  // get derivatives of base functions for columns
1094  auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1095 
1096  // iterate column base functions
1097  for (int cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1098 
1099  // integrate block local stiffens matrix
1100  t_m(i, k) += t_row(i, k, l) * t_col_diff_base(l);
1101 
1102  // move to next column base function
1103  ++t_col_diff_base;
1104 
1105  // move to next block of local stiffens matrix
1106  ++t_m;
1107  }
1108 
1109  // move to next row base function
1110  ++t_row_diff_base;
1111  }
1112 
1113  for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1114  ++t_row_diff_base;
1115 
1116  // move to next integration weight
1117  ++t_w;
1118  ++t_D;
1119  ++t_coords;
1120  }
1121  }
1122 
1124 }
1125 
1126 template <int SPACE_DIM, typename OpBase>
1128  EntitiesFieldData::EntData &row_data,
1129  EntitiesFieldData::EntData &col_data) {
1131 
1132  auto t_w = this->getFTensor0IntegrationWeight();
1133 
1134  size_t nb_base_functions = row_data.getN().size2() / 3;
1135  auto t_row_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
1136 
1137  const double alpha_constant = alphaConstant();
1138 
1139  for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1140  const double alpha = alpha_constant * this->getMeasure() * t_w;
1141 
1142  size_t rr = 0;
1143  for (; rr != OpBase::nbRows; ++rr) {
1144  const double t_row_div_base = t_row_diff_base(i, i);
1145  auto t_col_base = col_data.getFTensor0N(gg, 0);
1146  for (size_t cc = 0; cc != OpBase::nbCols; ++cc) {
1147  this->locMat(rr, cc) += alpha * t_row_div_base * t_col_base;
1148  ++t_col_base;
1149  }
1150  ++t_row_diff_base;
1151  }
1152  for (; rr < nb_base_functions; ++rr)
1153  ++t_row_diff_base;
1154 
1155  ++t_w;
1156  }
1157 
1159 }
1160 
1161 template <int SPACE_DIM, typename OpBase, CoordinateTypes CoordSys>
1164  EntitiesFieldData::EntData &row_data,
1165  EntitiesFieldData::EntData &col_data) {
1167 
1168  auto t_w = this->getFTensor0IntegrationWeight();
1169  auto t_coords = this->getFTensor1CoordsAtGaussPts();
1170 
1171  size_t nb_base_functions = row_data.getN().size2() / 3;
1172  auto t_row_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
1173  auto t_row_base = row_data.getFTensor1N<3>();
1174  const double alpha_constant = alphaConstant() * this->getMeasure();
1175  for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1176 
1177  const double alpha =
1178  alpha_constant * t_w * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1179 
1180  size_t rr = 0;
1181  for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1182  auto t_mat_diag = getFTensor1FromArrayDiag<SPACE_DIM, SPACE_DIM>(
1183  this->locMat, SPACE_DIM * rr);
1184  const double t_row_div_base = t_row_diff_base(i, i);
1185  auto t_col_base = col_data.getFTensor0N(gg, 0);
1186 
1187  for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1188  t_mat_diag(i) += alpha * t_row_div_base * t_col_base;
1189  if constexpr (CoordSys == CYLINDRICAL) {
1190  t_mat_diag(i) += t_row_base(0) * (alpha / t_coords(0)) * t_col_base;
1191  }
1192  ++t_col_base;
1193  ++t_mat_diag;
1194  }
1195  ++t_row_base;
1196  ++t_row_diff_base;
1197  }
1198  for (; rr < nb_base_functions; ++rr) {
1199  ++t_row_diff_base;
1200  ++t_row_base;
1201  }
1202 
1203  ++t_w;
1204  ++t_coords;
1205  }
1206 
1208 }
1209 
1210 template <int SPACE_DIM, typename OpBase, CoordinateTypes COORDINATE_SYSTEM>
1213  EntitiesFieldData::EntData &row_data,
1214  EntitiesFieldData::EntData &col_data) {
1216 
1217 #ifndef NDEBUG
1218  if (OpBase::locMat.size2() % SPACE_DIM)
1219  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1220  "Number of rows in matrix should be multiple of space dimensions");
1221 #endif
1222 
1223  // When we move to C++17 add if constexpr()
1224  if constexpr (COORDINATE_SYSTEM == POLAR || COORDINATE_SYSTEM == SPHERICAL)
1225  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1226  "%s coordiante not implemented",
1227  CoordinateTypesNames[COORDINATE_SYSTEM]);
1228 
1229  auto t_w = this->getFTensor0IntegrationWeight();
1230  auto t_coords = this->getFTensor1CoordsAtGaussPts();
1231  size_t nb_base_functions_row = row_data.getN().size2();
1232  auto t_row_base = row_data.getFTensor0N();
1233  const double vol = this->getMeasure();
1234  for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1235 
1236  const double alpha =
1237  alphaConstant(t_coords(0), t_coords(1), t_coords(2)) * t_w * vol;
1238 
1239  size_t rr = 0;
1240  auto t_m = getFTensor1FromPtr<SPACE_DIM>(OpBase::locMat.data().data());
1241 
1242  // When we move to C++17 add if constexpr()
1243  if constexpr (COORDINATE_SYSTEM == CARTESIAN) {
1244  for (; rr != OpBase::nbRows; ++rr) {
1245  const double r_val = alpha * t_row_base;
1246  auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1247  for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1248  t_m(i) += r_val * t_col_diff_base(i);
1249  ++t_col_diff_base;
1250  ++t_m;
1251  }
1252  ++t_row_base;
1253  }
1254  }
1255 
1256  // When we move to C++17 add if constexpr()
1257  if constexpr (COORDINATE_SYSTEM == CYLINDRICAL) {
1258  for (; rr != OpBase::nbRows; ++rr) {
1259  const double r_val = alpha * t_row_base;
1260  auto t_col_base = col_data.getFTensor0N(gg, 0);
1261  auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1262  for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1263  t_m(i) += r_val * t_col_diff_base(i);
1264  t_m(0) += (r_val / t_coords(0)) * t_col_base;
1265  ++t_col_base;
1266  ++t_col_diff_base;
1267  ++t_m;
1268  }
1269  ++t_row_base;
1270  }
1271  }
1272 
1273  for (; rr < nb_base_functions_row; ++rr)
1274  ++t_row_base;
1275 
1276  ++t_w;
1277  ++t_coords;
1278  }
1279 
1281 }
1282 
1283 template <int SPACE_DIM, typename OpBase>
1286  EntitiesFieldData::EntData &row_data,
1287  EntitiesFieldData::EntData &col_data) {
1289 
1290  auto t_w = this->getFTensor0IntegrationWeight();
1291 
1292  size_t nb_base_functions = row_data.getN().size2() / 3;
1293  auto t_row_base = row_data.getFTensor1N<3>();
1294  auto &mat = this->locMat;
1295  const double alpha_constant = alphaConstant();
1296  for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1297 
1298  const double alpha = alpha_constant * this->getMeasure() * t_w;
1299 
1300  size_t rr = 0;
1301  for (; rr != OpBase::nbRows; ++rr) {
1302  auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1303  for (size_t cc = 0; cc != OpBase::nbCols; ++cc) {
1304  mat(rr, cc) += alpha * t_row_base(i) * t_col_diff_base(i);
1305  ++t_col_diff_base;
1306  }
1307  ++t_row_base;
1308  }
1309  for (; rr < nb_base_functions; ++rr)
1310  ++t_row_base;
1311 
1312  ++t_w;
1313  }
1314 
1316 }
1317 
1318 template <int SPACE_DIM, typename OpBase>
1321  EntitiesFieldData::EntData &row_data,
1322  EntitiesFieldData::EntData &col_data) {
1324 
1325  auto t_w = this->getFTensor0IntegrationWeight();
1326 
1327  size_t nb_base_functions = row_data.getN().size2();
1328  auto t_row_base = row_data.getFTensor0N();
1329 
1330  auto get_t_vec = [&](const int rr) {
1331  std::array<double *, SPACE_DIM> ptrs;
1332  for (auto i = 0; i != SPACE_DIM; ++i)
1333  ptrs[i] = &OpBase::locMat(rr + i, 0);
1335  };
1336 
1337  const double alpha_constant = alphaConstant();
1338  for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1339 
1340  const double alpha = alpha_constant * this->getMeasure() * t_w;
1341 
1342  size_t rr = 0;
1343  for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1344  auto t_vec = get_t_vec(SPACE_DIM * rr);
1345  auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1346  for (size_t cc = 0; cc != OpBase::nbCols; ++cc) {
1347  t_vec(i) += alpha * t_row_base * t_col_diff_base(i);
1348  ++t_col_diff_base;
1349  ++t_vec;
1350  }
1351  ++t_row_base;
1352  }
1353  for (; rr < nb_base_functions; ++rr)
1354  ++t_row_base;
1355 
1356  ++t_w;
1357  }
1358 
1360 }
1361 
1362 template <int SPACE_DIM, typename OpBase>
1364  EntitiesFieldData::EntData &row_data,
1365  EntitiesFieldData::EntData &col_data) {
1367 
1368  auto t_w = this->getFTensor0IntegrationWeight();
1369  auto t_coords = this->getFTensor1CoordsAtGaussPts();
1370 
1371  size_t nb_base_functions = row_data.getN().size2() / 3;
1372  auto t_row_base = row_data.getFTensor1N<3>();
1373  const double alpha_constant = alphaConstant() * this->getMeasure();
1374  for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1375 
1376  const double alpha =
1377  alpha_constant * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * t_w;
1378 
1379  size_t rr = 0;
1380  for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1381  auto t_mat_diag = getFTensor1FromArrayDiag<SPACE_DIM, SPACE_DIM>(
1382  this->locMat, SPACE_DIM * rr);
1383  auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1384 
1385  for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1386  t_mat_diag(i) += alpha * t_row_base(j) * t_col_diff_base(j);
1387  ++t_col_diff_base;
1388  ++t_mat_diag;
1389  }
1390 
1391  ++t_row_base;
1392  }
1393  for (; rr < nb_base_functions; ++rr)
1394  ++t_row_base;
1395 
1396  ++t_w;
1397  ++t_coords;
1398  }
1399 
1401 }
1402 
1403 } // namespace MoFEM
1404 
1405 #endif //__BILINEAR_FORMS_INTEGRATORS_IMPL_HPP__
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::OpMixTensorTimesGradImpl< SPACE_DIM, GAUSS, OpBase >::j
FTensor::Index< 'j', SPACE_DIM > j
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:414
MoFEM::OpMixTensorTimesGradImpl< SPACE_DIM, GAUSS, OpBase >::OpMixTensorTimesGradImpl
OpMixTensorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, ScalarFun beta_coeff, const bool assemble_transpose, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:401
MoFEM::OpBaseImpl::nbIntegrationPts
int nbIntegrationPts
number of integration points
Definition: FormsIntegrators.hpp:238
FTensor::Christof
Definition: Christof_value.hpp:9
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:128
MoFEM::OpMassCacheImpl< 1, FIELD_DIM, I, OpBase >::scalarBeta
double scalarBeta
Definition: BiLinearFormsIntegratorsImpl.hpp:177
MoFEM::EntitiesFieldData::EntData::getFieldEntities
const VectorFieldEntities & getFieldEntities() const
get field entities
Definition: EntitiesFieldData.hpp:1267
MoFEM::OpBaseImpl::locMat
MatrixDouble locMat
local entity block matrix
Definition: FormsIntegrators.hpp:249
MoFEM::OpMassCacheImpl< 1, FIELD_DIM, I, OpBase >::OpMassCacheImpl
OpMassCacheImpl(CacheMatsTypeType cache, const std::string row_field_name, const std::string col_field_name, const double beta, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:165
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
CARTESIAN
@ CARTESIAN
Definition: definitions.h:128
MoFEM::OpMassImpl< 1, 1, GAUSS, OpBase >::OpMassImpl
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr, boost::shared_ptr< MatrixDouble > cache_mat=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:69
MoFEM::OpGradSymTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::OpGradSymTensorGradImpl
OpGradSymTensorGradImpl(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< MatrixDouble > mat_D, boost::shared_ptr< Range > ents_ptr=nullptr, ScalarFun beta=scalar_fun_one)
Definition: BiLinearFormsIntegratorsImpl.hpp:189
MoFEM::OpGradGradImpl< 1, 1, SPACE_DIM, GAUSS, OpBase >::OpGradGradImpl
OpGradGradImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:28
MoFEM::OpMassImpl< 3, FIELD_DIM, GAUSS, OpBase >::OpMassImpl
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:110
MoFEM::OpGradGradSymTensorGradGradImpl< 1, 1, SPACE_DIM, S, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:236
MoFEM::OpMassImpl< 3, 9, GAUSS, OpBase >::OpMassImpl
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:144
MoFEM::OpMixDivTimesVecImpl< SPACE_DIM, GAUSS, OpBase, CoordSys >::OpMixDivTimesVecImpl
OpMixDivTimesVecImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:284
MoFEM::OpMixDivTimesVecImpl< SPACE_DIM, GAUSS, OpBase, CoordSys >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:305
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::OpMixDivTimesScalarImpl< SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:270
MoFEM::ConstantFun
boost::function< double()> ConstantFun
Constant function type.
Definition: FormsIntegrators.hpp:166
FTensor::Tensor2_symmetric
Definition: Tensor2_symmetric_value.hpp:13
OpBase
OpBaseImpl< PETSC, EdgeEleOp > OpBase
Definition: radiation.cpp:29
BASE_DIM
constexpr int BASE_DIM
Definition: dg_projection.cpp:15
CoordinateTypes
CoordinateTypes
Coodinate system.
Definition: definitions.h:127
MoFEM::OpGradSymTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::matD
boost::shared_ptr< MatrixDouble > matD
Definition: BiLinearFormsIntegratorsImpl.hpp:201
FIELD_DIM
constexpr int FIELD_DIM
Definition: child_and_parent.cpp:15
MoFEM::OpGradSymTensorGradImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:183
I
constexpr IntegrationType I
Definition: operators_tests.cpp:31
MoFEM::scalar_fun_one
double scalar_fun_one(const double, const double, const double)
Definition: FormsIntegrators.hpp:146
MoFEM::OpMassCacheImpl< 1, FIELD_DIM, I, OpBase >::cacheLocMats
CacheMatsTypeType cacheLocMats
Definition: BiLinearFormsIntegratorsImpl.hpp:178
MoFEM::EntitiesFieldData::EntData::getFTensor0N
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
Definition: EntitiesFieldData.hpp:1492
MoFEM::OpMixScalarTimesDivImpl< SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM >::OpMixScalarTimesDivImpl
OpMixScalarTimesDivImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun alpha_fun=[](double, double, double) constexpr { return 1;}, const bool assemble_transpose=false, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:321
MoFEM::OpGradSymTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:188
MoFEM::OpGradGradImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:23
MoFEM::EntitiesFieldData::EntData::getFTensor2DiffN
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN(FieldApproximationBase base)
Get derivatives of base functions for Hdiv space.
Definition: EntitiesFieldData.cpp:660
MoFEM::OpBaseImpl
Definition: FormsIntegrators.hpp:178
MoFEM::OpGradGradImpl< 1, 1, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:27
CoordinateTypesNames
const static char *const CoordinateTypesNames[]
Coordinate system names.
Definition: definitions.h:121
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
FTensor::Tensor3
Definition: Tensor3_value.hpp:12
MoFEM::OpMixVectorTimesGradImpl< 3, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::alphaConstant
ConstantFun alphaConstant
Definition: BiLinearFormsIntegratorsImpl.hpp:359
MoFEM::OpGradGradSymTensorGradGradImpl< 1, 1, SPACE_DIM, S, GAUSS, OpBase >::OpGradGradSymTensorGradGradImpl
OpGradGradSymTensorGradGradImpl(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< MatrixDouble > mat_D, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:237
MoFEM::OpMassImpl< 1, FIELD_DIM, GAUSS, OpBase >::iNtegrate
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Integrate grad-grad operator.
Definition: BiLinearFormsIntegratorsImpl.hpp:101
MoFEM::OpGradTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::OpGradTensorGradImpl
OpGradTensorGradImpl(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< MatrixDouble > mat_D, ScalarFun beta=scalar_fun_one)
Definition: BiLinearFormsIntegratorsImpl.hpp:215
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
SPACE_DIM
constexpr int SPACE_DIM
Definition: child_and_parent.cpp:16
MoFEM::OpMassImpl< 3, 4, GAUSS, OpBase >::OpMassImpl
OpMassImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:127
a
constexpr double a
Definition: approx_sphere.cpp:30
POLAR
@ POLAR
Definition: definitions.h:129
MoFEM::OpMixVectorTimesGradImpl< 3, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:358
MoFEM::OpMixTensorTimesGradImpl< SPACE_DIM, GAUSS, OpBase >::OpMixTensorTimesGradImpl
OpMixTensorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:391
MoFEM::OpGradTensorGradImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:209
double
MoFEM::OpMassImpl< 3, FIELD_DIM, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition: BiLinearFormsIntegratorsImpl.hpp:120
MoFEM::OpMassImpl< 3, 9, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition: BiLinearFormsIntegratorsImpl.hpp:154
MoFEM::OpMixDivTimesVecImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:278
MoFEM::OpGradSymTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition: BiLinearFormsIntegratorsImpl.hpp:202
MoFEM::EntitiesFieldData::EntData::getIndices
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Definition: EntitiesFieldData.hpp:1204
MoFEM::OpGradTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition: BiLinearFormsIntegratorsImpl.hpp:224
MoFEM::OpMixVectorTimesGradImpl< 3, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::OpMixVectorTimesGradImpl
OpMixVectorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose=false, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:346
MoFEM::OpGradGradImpl< 1, FIELD_DIM, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:46
SPHERICAL
@ SPHERICAL
Definition: definitions.h:131
MoFEM::GAUSS
@ GAUSS
Definition: FormsIntegrators.hpp:136
MoFEM::OpBaseImpl::colSide
int colSide
column side number
Definition: FormsIntegrators.hpp:242
MoFEM::OpGradGradImpl< 1, FIELD_DIM, SPACE_DIM, GAUSS, OpBase >::OpGradGradImpl
OpGradGradImpl(const std::string row_field_name, const std::string col_field_name, ScalarFun beta=scalar_fun_one, boost::shared_ptr< Range > ents_ptr=nullptr)
Definition: BiLinearFormsIntegratorsImpl.hpp:47
MoFEM::OpGradGradImpl< 1, 1, SPACE_DIM, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition: BiLinearFormsIntegratorsImpl.hpp:39
MoFEM::ScalarFun
boost::function< double(const double, const double, const double)> ScalarFun
Scalar function type.
Definition: FormsIntegrators.hpp:144
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MoFEM::OpGradTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:214
MoFEM::OpMixScalarTimesDivImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:316
EntData
EntitiesFieldData::EntData EntData
Definition: child_and_parent.cpp:37
FTensor::Index< 'i', SPACE_DIM >
MoFEM::IntegrationType
IntegrationType
Form integrator integration types.
Definition: FormsIntegrators.hpp:136
MoFEM::OpBaseImpl::nbRowBaseFunctions
int nbRowBaseFunctions
number or row base functions
Definition: FormsIntegrators.hpp:239
MoFEM::OpBaseImpl::rowSide
int rowSide
row side number
Definition: FormsIntegrators.hpp:241
MoFEM::OpBaseImpl::nbCols
int nbCols
number if dof on column
Definition: FormsIntegrators.hpp:237
MoFEM::OpMixVectorTimesGradImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:341
MoFEM::OpMassImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:64
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::EntitiesFieldData::EntData::getFTensor1DiffN
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
Definition: EntitiesFieldData.cpp:526
MoFEM::OpBaseImpl::nbRows
int nbRows
number of dofs on rows
Definition: FormsIntegrators.hpp:236
CYLINDRICAL
@ CYLINDRICAL
Definition: definitions.h:130
MoFEM::OpMixDivTimesScalarImpl< SPACE_DIM, GAUSS, OpBase >::alphaConstant
ConstantFun alphaConstant
Definition: BiLinearFormsIntegratorsImpl.hpp:271
MoFEM::OpMixTensorTimesGradImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:386
MoFEM::OpMassImpl< 1, 1, GAUSS, OpBase >::iNtegrate
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Integrate grad-grad operator.
Definition: BiLinearFormsIntegratorsImpl.hpp:85
MoFEM::OpMixVectorTimesGradImpl< 1, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::OpMixVectorTimesGradImpl
OpMixVectorTimesGradImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose=false, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:367
MoFEM::EntitiesFieldData::EntData::getN
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
Definition: EntitiesFieldData.hpp:1308
MoFEM::OpMixVectorTimesGradImpl< 1, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::alphaConstant
ConstantFun alphaConstant
Definition: BiLinearFormsIntegratorsImpl.hpp:380
MoFEM::OpMassCacheImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:160
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
MoFEM::EntitiesFieldData::EntData::getFTensor1N
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
Definition: EntitiesFieldData.cpp:640
MoFEM::OpMixVectorTimesGradImpl< 1, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:379
MoFEM::OpMixDivTimesScalarImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:254
MoFEM::OpGradGradSymTensorGradGradImpl
Definition: BiLinearFormsIntegratorsImpl.hpp:231
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::OpGradGradImpl< 1, FIELD_DIM, SPACE_DIM, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition: BiLinearFormsIntegratorsImpl.hpp:58
MoFEM::OpMixDivTimesScalarImpl< SPACE_DIM, GAUSS, OpBase >::OpMixDivTimesScalarImpl
OpMixDivTimesScalarImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, const bool assemble_transpose=false, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:258
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
MoFEM::OpGradGradSymTensorGradGradImpl< 1, 1, SPACE_DIM, S, GAUSS, OpBase >::matD
boost::shared_ptr< MatrixDouble > matD
Definition: BiLinearFormsIntegratorsImpl.hpp:248
MoFEM::OpMixScalarTimesDivImpl< SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM >::alphaConstant
ScalarFun alphaConstant
Definition: BiLinearFormsIntegratorsImpl.hpp:334
MoFEM::CacheMatsTypeType
std::map< std::pair< int, int >, boost::shared_ptr< MatrixDouble > > CacheMatsTypeType
Definition: BiLinearFormsIntegratorsImpl.hpp:19
MoFEM::OpBaseImpl::rowType
EntityType rowType
row type
Definition: FormsIntegrators.hpp:243
MoFEM::OpMixTensorTimesGradImpl< SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:413
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEM::OpMixDivTimesVecImpl< SPACE_DIM, GAUSS, OpBase, CoordSys >::OpMixDivTimesVecImpl
OpMixDivTimesVecImpl(const std::string row_field_name, const std::string col_field_name, ConstantFun alpha_fun, ScalarFun beta_fun, const bool assemble_transpose, const bool only_transpose=false)
Definition: BiLinearFormsIntegratorsImpl.hpp:294
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEM::OpGradTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::matD
boost::shared_ptr< MatrixDouble > matD
Definition: BiLinearFormsIntegratorsImpl.hpp:223
MoFEM::OpMassImpl< 3, 4, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition: BiLinearFormsIntegratorsImpl.hpp:137
MoFEM::OpMixScalarTimesDivImpl< SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition: BiLinearFormsIntegratorsImpl.hpp:333
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::OpBaseImpl::colType
EntityType colType
column type
Definition: FormsIntegrators.hpp:244
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21