v0.16.3
Loading...
Searching...
No Matches
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
16namespace MoFEM {
17
18using CacheMatsTypeType =
19 std::map<std::pair<int, int>, boost::shared_ptr<MatrixDouble>>;
20
21template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
22 typename OpBase>
23struct OpGradGradImpl {};
24
25template <int SPACE_DIM, typename OpBase>
26struct 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,
30 ScalarFun beta = scalar_fun_one,
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
38protected:
39 ScalarFun betaCoeff;
40 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
41 EntitiesFieldData::EntData &col_data);
42};
43
44template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
45struct OpGradGradImpl<1, FIELD_DIM, SPACE_DIM, GAUSS, OpBase> : public 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,
49 ScalarFun beta = scalar_fun_one,
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
57protected:
58 ScalarFun betaCoeff;
59 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
60 EntitiesFieldData::EntData &col_data);
61};
62
63template <int BASE_DIM, int FIELD_DIM, IntegrationType I, typename OpBase>
64struct OpMassImpl {};
65
66template <typename OpBase>
67struct OpMassImpl<1, 1, GAUSS, OpBase> : public OpBase {
68
69 OpMassImpl(
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
80protected:
81 using OpBase::OpBase;
82
83 ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
84 MoFEMErrorCode integrateImpl(EntitiesFieldData::EntData &row_data,
85 EntitiesFieldData::EntData &col_data,
86 double vol);
87 inline MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
88 EntitiesFieldData::EntData &col_data) {
89 return integrateImpl(row_data, col_data, OpBase::getMeasure());
90 }
91};
92
93template <int FIELD_DIM, typename OpBase>
94struct OpMassImpl<1, FIELD_DIM, GAUSS, OpBase>
95 : public OpMassImpl<1, 1, GAUSS, OpBase> {
96 using OpMassImpl<1, 1, GAUSS, OpBase>::OpMassImpl;
97
98protected:
99 using OpMassImpl<1, 1, GAUSS, OpBase>::betaCoeff;
100 MoFEMErrorCode integrateImpl(EntitiesFieldData::EntData &row_data,
101 EntitiesFieldData::EntData &col_data,
102 double vol);
103 inline MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
104 EntitiesFieldData::EntData &col_data) {
105 return integrateImpl(row_data, col_data, OpBase::getMeasure());
106 }
110};
111
112template <int FIELD_DIM, typename OpBase>
113struct OpMassImpl<3, FIELD_DIM, GAUSS, OpBase> : public OpBase {
114
115 OpMassImpl(const std::string row_field_name, const std::string col_field_name,
116 ScalarFun beta = scalar_fun_one,
117 boost::shared_ptr<Range> ents_ptr = nullptr)
118 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
119 betaCoeff(beta) {
120 if (row_field_name == col_field_name)
121 this->sYmm = true;
122 }
123
124protected:
125 ScalarFun betaCoeff;
126 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
127 EntitiesFieldData::EntData &col_data);
128};
129
130template <typename OpBase>
131struct OpMassImpl<3, 4, GAUSS, OpBase> : public OpBase {
132 OpMassImpl(const std::string row_field_name, const std::string col_field_name,
133 ScalarFun beta = scalar_fun_one,
134 boost::shared_ptr<Range> ents_ptr = nullptr)
135 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
136 betaCoeff(beta) {
137 if (row_field_name == col_field_name)
138 this->sYmm = true;
139 }
140
141protected:
142 ScalarFun betaCoeff;
143 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
144 EntitiesFieldData::EntData &col_data);
145};
146
147template <typename OpBase>
148struct OpMassImpl<3, 9, GAUSS, OpBase> : public OpBase {
149 OpMassImpl(const std::string row_field_name, const std::string col_field_name,
150 ScalarFun beta = scalar_fun_one,
151 boost::shared_ptr<Range> ents_ptr = nullptr)
152 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
153 betaCoeff(beta) {
154 if (row_field_name == col_field_name)
155 this->sYmm = true;
156 }
157
158protected:
159 ScalarFun betaCoeff;
160 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
161 EntitiesFieldData::EntData &col_data);
162};
163
164template <typename OpBase>
165struct OpMassImpl<9, 9, GAUSS, OpBase> : public OpBase {
166 OpMassImpl(const std::string row_field_name, const std::string col_field_name,
167 ScalarFun beta = scalar_fun_one,
168 boost::shared_ptr<Range> ents_ptr = nullptr)
169 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
170 betaCoeff(beta) {
171 if (row_field_name == col_field_name)
172 this->sYmm = true;
173 }
174
175protected:
176 ScalarFun betaCoeff;
177 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
178 EntitiesFieldData::EntData &col_data);
179};
180
181template <int BASE_DIM, int FIELD_DIM, IntegrationType I, typename OpBase>
182struct OpMassCacheImpl;
183
184template <int FIELD_DIM, IntegrationType I, typename OpBase>
185struct OpMassCacheImpl<1, FIELD_DIM, I, OpBase>
186 : public OpMassImpl<1, FIELD_DIM, I, OpBase> {
187 OpMassCacheImpl(CacheMatsTypeType cache, const std::string row_field_name,
188 const std::string col_field_name, const double beta,
189 boost::shared_ptr<Range> ents_ptr = nullptr)
190 : OpMassImpl<1, FIELD_DIM, I, OpBase>(
191 row_field_name, col_field_name,
192 [](double, double, double) constexpr { return 1; }, ents_ptr),
193 cacheLocMats(cache), scalarBeta(beta) {}
194
195 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
196 EntitiesFieldData::EntData &col_data);
197
198protected:
199 double scalarBeta;
200 CacheMatsTypeType cacheLocMats;
201};
202
203template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
204 typename OpBase>
205struct OpGradSymTensorGradImpl {};
206
207template <int SPACE_DIM, int S, typename OpBase>
208struct OpGradSymTensorGradImpl<1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase>
209 : public OpBase {
210 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
211 OpGradSymTensorGradImpl(const std::string row_field_name,
212 const std::string col_field_name,
213 boost::shared_ptr<MatrixDouble> mat_D,
214 boost::shared_ptr<Range> ents_ptr = nullptr,
215 ScalarFun beta = scalar_fun_one)
216 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
217 matD(mat_D), betaCoeff(beta) {
218 if (row_field_name == col_field_name)
219 this->sYmm = true;
220 }
221
222protected:
223 boost::shared_ptr<MatrixDouble> matD;
224 ScalarFun betaCoeff;
225 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
226 EntitiesFieldData::EntData &col_data);
227};
228
229template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
230 typename OpBase>
231struct OpGradTensorGradImpl {};
232
233template <int SPACE_DIM, int S, typename OpBase>
234struct OpGradTensorGradImpl<1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase>
235 : public OpBase {
236 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
237 OpGradTensorGradImpl(const std::string row_field_name,
238 const std::string col_field_name,
239 boost::shared_ptr<MatrixDouble> mat_D,
240 ScalarFun beta = scalar_fun_one)
241 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL), matD(mat_D),
242 betaCoeff(beta) {}
243
244protected:
245 boost::shared_ptr<MatrixDouble> matD;
246 ScalarFun betaCoeff;
247 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
248 EntitiesFieldData::EntData &col_data);
249};
250
251template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
252 typename OpBase>
253struct OpGradGradSymTensorGradGradImpl {};
254
255template <int FIELD_DIM, int SPACE_DIM, int S, typename OpBase>
256struct OpGradGradSymTensorGradGradImpl<1, FIELD_DIM, SPACE_DIM, S, GAUSS,
257 OpBase> : public OpBase {
258 OpGradGradSymTensorGradGradImpl(const std::string row_field_name,
259 const std::string col_field_name,
260 boost::shared_ptr<MatrixDouble> mat_D,
261 boost::shared_ptr<Range> ents_ptr = nullptr)
262 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
263 matD(mat_D) {
264 if (row_field_name == col_field_name)
265 this->sYmm = true;
266 }
267
268protected:
269 boost::shared_ptr<MatrixDouble> matD;
270 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
271 EntitiesFieldData::EntData &col_data);
272};
273
274template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
275 typename OpBase>
276struct OpGradGradTensorGradGradImpl {};
277
278template <int FIELD_DIM, int SPACE_DIM, int S, typename OpBase>
279struct OpGradGradTensorGradGradImpl<1, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase>
280 : public OpBase {
281 OpGradGradTensorGradGradImpl(const std::string row_field_name,
282 const std::string col_field_name,
283 boost::shared_ptr<MatrixDouble> mat_D,
284 boost::shared_ptr<Range> ents_ptr = nullptr)
285 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL, ents_ptr),
286 matD(mat_D) {
287 if (row_field_name == col_field_name)
288 this->sYmm = true;
289 }
290
291protected:
292 boost::shared_ptr<MatrixDouble> matD;
293 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
294 EntitiesFieldData::EntData &col_data);
295};
296
297template <int SPACE_DIM, IntegrationType I, typename OpBase>
298struct OpMixDivTimesScalarImpl {};
299
300template <int SPACE_DIM, typename OpBase>
301struct OpMixDivTimesScalarImpl<SPACE_DIM, GAUSS, OpBase> : public OpBase {
302 OpMixDivTimesScalarImpl(const std::string row_field_name,
303 const std::string col_field_name,
304 ConstantFun alpha_fun,
305 const bool assemble_transpose = false,
306 const bool only_transpose = false)
307 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
308 alphaConstant(alpha_fun) {
309 this->assembleTranspose = assemble_transpose;
310 this->onlyTranspose = only_transpose;
311 }
312
313protected:
314 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
315 ConstantFun alphaConstant;
316 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
317 EntitiesFieldData::EntData &col_data);
318};
319
320template <int SPACE_DIM, IntegrationType I, typename OpBase,
321 CoordinateTypes CoordSys>
322struct OpMixDivTimesVecImpl {};
323
324template <int SPACE_DIM, typename OpBase, CoordinateTypes CoordSys>
325struct OpMixDivTimesVecImpl<SPACE_DIM, GAUSS, OpBase, CoordSys>
326 : public OpBase {
327
328 OpMixDivTimesVecImpl(const std::string row_field_name,
329 const std::string col_field_name, ConstantFun alpha_fun,
330 const bool assemble_transpose,
331 const bool only_transpose = false)
332 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
333 alphaConstant(alpha_fun) {
334 this->assembleTranspose = assemble_transpose;
335 this->onlyTranspose = only_transpose;
336 }
337
338 OpMixDivTimesVecImpl(const std::string row_field_name,
339 const std::string col_field_name, ConstantFun alpha_fun,
340 ScalarFun beta_fun, const bool assemble_transpose,
341 const bool only_transpose = false)
342 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
343 alphaConstant(alpha_fun), betaCoeff(beta_fun) {
344 this->assembleTranspose = assemble_transpose;
345 this->onlyTranspose = only_transpose;
346 }
347
348protected:
349 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
350
351 ConstantFun alphaConstant = []() constexpr { return 1; };
352 ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
353
354 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
355 EntitiesFieldData::EntData &col_data);
356};
357
358template <int SPACE_DIM, IntegrationType I, typename OpBase,
359 CoordinateTypes COORDINATE_SYSTEM>
360struct OpMixScalarTimesDivImpl {};
361
362template <int SPACE_DIM, typename OpBase, CoordinateTypes COORDINATE_SYSTEM>
363struct OpMixScalarTimesDivImpl<SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM>
364 : public OpBase {
365 OpMixScalarTimesDivImpl(
366 const std::string row_field_name, const std::string col_field_name,
367 ScalarFun alpha_fun = [](double, double, double) constexpr { return 1; },
368 const bool assemble_transpose = false, const bool only_transpose = false)
369 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
370 alphaConstant(alpha_fun) {
371 this->assembleTranspose = assemble_transpose;
372 this->onlyTranspose = only_transpose;
373 this->sYmm = false;
374 }
375
376protected:
377 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
378 ScalarFun alphaConstant;
379 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
380 EntitiesFieldData::EntData &col_data);
381};
382
383template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
384 typename OpBase>
385struct OpMixVectorTimesGradImpl;
386
387template <int SPACE_DIM, typename OpBase>
388struct OpMixVectorTimesGradImpl<3, SPACE_DIM, SPACE_DIM, GAUSS, OpBase>
389 : public OpBase {
390 OpMixVectorTimesGradImpl(const std::string row_field_name,
391 const std::string col_field_name,
392 ConstantFun alpha_fun,
393 const bool assemble_transpose = false,
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
401protected:
402 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
403 ConstantFun alphaConstant;
404 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
405 EntitiesFieldData::EntData &col_data);
406};
407
408template <int SPACE_DIM, typename OpBase>
409struct OpMixVectorTimesGradImpl<1, SPACE_DIM, SPACE_DIM, GAUSS, OpBase>
410 : public OpBase {
411 OpMixVectorTimesGradImpl(const std::string row_field_name,
412 const std::string col_field_name,
413 ConstantFun alpha_fun,
414 const bool assemble_transpose = false,
415 const bool only_transpose = false)
416 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
417 alphaConstant(alpha_fun) {
418 this->assembleTranspose = assemble_transpose;
419 this->onlyTranspose = only_transpose;
420 }
421
422protected:
423 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
424 ConstantFun alphaConstant;
425 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
426 EntitiesFieldData::EntData &col_data);
427};
428
429template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
430 typename OpBase>
431struct OpBaseTensorBaseImpl {};
432
433template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, typename OpBase>
434struct OpBaseTensorBaseImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase>
435 : public OpBase {
436 static_assert(BASE_DIM == 2 || BASE_DIM == 3 || BASE_DIM == 4 || BASE_DIM == 9,
437 "Expected vector or square tensor column bases");
438 static_assert(FIELD_DIM > 0 && (SPACE_DIM == 2 || SPACE_DIM == 3),
439 "Expected a positive row field size and 2D or 3D tensors");
440
441 OpBaseTensorBaseImpl(
442 const std::string row_field_name, const std::string col_field_name,
443 boost::shared_ptr<MatrixDouble> mat_tensor,
444 ConstantFun alpha_fun = []() constexpr { return 1.; },
445 const bool assemble_transpose = false, const bool only_transpose = false)
446 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
447 matTensor(mat_tensor), alphaConstant(alpha_fun) {
448 this->sYmm = false;
449 this->assembleTranspose = assemble_transpose;
450 this->onlyTranspose = only_transpose;
451 }
452
453protected:
454 boost::shared_ptr<MatrixDouble> matTensor;
455 ConstantFun alphaConstant;
456 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
457 EntitiesFieldData::EntData &col_data) override;
458};
459
460template <int SPACE_DIM, IntegrationType I, typename OpBase>
461struct OpMixTensorTimesGradImpl {};
462
463template <int SPACE_DIM, typename OpBase>
464struct OpMixTensorTimesGradImpl<SPACE_DIM, GAUSS, OpBase> : public OpBase {
465
466 OpMixTensorTimesGradImpl(const std::string row_field_name,
467 const std::string col_field_name,
468 ConstantFun alpha_fun, const bool assemble_transpose,
469 const bool only_transpose = false)
470 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
471 alphaConstant(alpha_fun) {
472 this->assembleTranspose = assemble_transpose;
473 this->onlyTranspose = only_transpose;
474 }
475
476 OpMixTensorTimesGradImpl(const std::string row_field_name,
477 const std::string col_field_name,
478 ConstantFun alpha_fun, ScalarFun beta_coeff,
479 const bool assemble_transpose,
480 const bool only_transpose = false)
481 : OpBase(row_field_name, col_field_name, OpBase::OPROWCOL),
482 alphaConstant(alpha_fun), betaCoeff(beta_coeff) {
483 this->assembleTranspose = assemble_transpose;
484 this->onlyTranspose = only_transpose;
485 }
486
487protected:
488 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
489 FTensor::Index<'j', SPACE_DIM> j; ///< summit Index
490 ConstantFun alphaConstant = []() constexpr { return 1; };
491 ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
492 MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data,
493 EntitiesFieldData::EntData &col_data);
494};
495
496template <int SPACE_DIM, typename OpBase>
497MoFEMErrorCode OpGradGradImpl<1, 1, SPACE_DIM, GAUSS, OpBase>::iNtegrate(
498 EntitiesFieldData::EntData &row_data,
499 EntitiesFieldData::EntData &col_data) {
500 MoFEMFunctionBegin;
501
502 // get element volume
503 const double vol = OpBase::getMeasure();
504 // get integration weights
505 auto t_w = OpBase::getFTensor0IntegrationWeight();
506 // get base function gradient on rows
507 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
508 // get coordinate at integration points
509 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
510
511 // loop over integration points
512 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
513 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
514 // take into account Jacobian
515 const double alpha = t_w * beta;
516 // loop over ros base functions
517 int rr = 0;
518 for (; rr != OpBase::nbRows; rr++) {
519 // get column base functions gradient at gauss point gg
520 auto t_col_grad = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
521 // loop over columns
522 for (int cc = 0; cc != OpBase::nbCols; cc++) {
523 // calculate element of local matrix
524 OpBase::locMat(rr, cc) += alpha * (t_row_grad(i) * t_col_grad(i));
525 ++t_col_grad; // move to another gradient of base function
526 // on column
527 }
528 ++t_row_grad; // move to another element of gradient of base
529 // function on row
530 }
531 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
532 ++t_row_grad;
533
534 ++t_coords;
535 ++t_w; // move to another integration weight
536 }
537 MoFEMFunctionReturn(0);
538}
539
540template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
541MoFEMErrorCode
542OpGradGradImpl<1, FIELD_DIM, SPACE_DIM, GAUSS, OpBase>::iNtegrate(
543 EntitiesFieldData::EntData &row_data,
544 EntitiesFieldData::EntData &col_data) {
545 MoFEMFunctionBegin;
546
547 // get element volume
548 const double vol = OpBase::getMeasure();
549 // get integration weights
550 auto t_w = OpBase::getFTensor0IntegrationWeight();
551 // get base function gradient on rows
552 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
553 // get coordinate at integration points
554 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
555
556 FTensor::Index<'i', SPACE_DIM> i;
557 FTensor::Index<'j', FIELD_DIM> j;
558
559 auto get_t_vec = [&](const int rr) {
560 std::array<double *, FIELD_DIM> ptrs;
561 for (auto i = 0; i != FIELD_DIM; ++i)
562 ptrs[i] = &OpBase::locMat(rr + i, i);
563 return FTensor::Tensor1<FTensor::PackPtr<double *, FIELD_DIM>, FIELD_DIM>(
564 ptrs);
565 };
566
567 // loop over integration points
568 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
569 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
570 // take into account Jacobian
571 const double alpha = t_w * beta;
572 // loop over ros base functions
573 int rr = 0;
574 for (; rr != OpBase::nbRows / FIELD_DIM; rr++) {
575 // get diag vec
576 auto t_vec = get_t_vec(rr * FIELD_DIM);
577 // get column base functions gradient at gauss point gg
578 auto t_col_grad = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
579 // loop over columns
580 for (int cc = 0; cc != OpBase::nbCols / FIELD_DIM; cc++) {
581 // calculate element of local matrix
582 t_vec(j) += alpha * (t_row_grad(i) * t_col_grad(i));
583 ++t_col_grad; // move to another gradient of base function
584 // on column
585 ++t_vec;
586 }
587 ++t_row_grad; // move to another element of gradient of base
588 // function on row
589 }
590 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
591 ++t_row_grad;
592
593 ++t_coords;
594 ++t_w; // move to another integration weight
595 }
596 MoFEMFunctionReturn(0);
597}
598
599template <typename OpBase>
600MoFEMErrorCode OpMassImpl<1, 1, GAUSS, OpBase>::integrateImpl(
601 EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data,
602 double vol) {
603 MoFEMFunctionBegin;
604
605#ifndef NDEBUG
606 auto log_error = [&]() {
607 MOFEM_LOG("SELF", Sev::error) << "Row side " << OpBase::rowSide << " "
608 << CN::EntityTypeName(OpBase::rowType);
609 MOFEM_LOG("SELF", Sev::error) << "Col side " << OpBase::colSide << " "
610 << CN::EntityTypeName(OpBase::colType);
611 };
612
613 if (row_data.getN().size2() < OpBase::nbRows) {
614 log_error();
615 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
616 "Wrong number of base functions on rows %zu < %d",
617 row_data.getN().size2(), OpBase::nbRows);
618 }
619 if (col_data.getN().size2() < OpBase::nbCols) {
620 log_error();
621 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
622 "Wrong number of base functions on cols %zu < %d",
623 col_data.getN().size2(), OpBase::nbCols);
624 }
625 if (row_data.getN().size1() != OpBase::nbIntegrationPts) {
626 log_error();
627 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
628 "Wrong number of integration points on rows %zu != %d",
629 row_data.getN().size1(), OpBase::nbIntegrationPts);
630 }
631 if (col_data.getN().size1() != OpBase::nbIntegrationPts) {
632 log_error();
633 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
634 "Wrong number of integration points on cols %zu != %d",
635 col_data.getN().size1(), OpBase::nbIntegrationPts);
636 }
637#endif
638
639 // get integration weights
640 auto t_w = OpBase::getFTensor0IntegrationWeight();
641 // get base function gradient on rows
642 auto t_row_base = row_data.getFTensor0N();
643 // get coordinate at integration points
644 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
645 // loop over integration points
646 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
647 const double beta = betaCoeff(t_coords(0), t_coords(1), t_coords(2));
648 // take into account Jacobian
649 const double alpha = t_w * beta;
650 // loop over rows base functions
651 auto a_mat_ptr = &*OpBase::locMat.data().begin();
652 int rr = 0;
653 for (; rr != OpBase::nbRows; rr++) {
654 // get column base functions gradient at gauss point gg
655 auto t_col_base = col_data.getFTensor0N(gg, 0);
656 // loop over columns
657 for (int cc = 0; cc != OpBase::nbCols; cc++) {
658 // calculate element of local matrix
659 *a_mat_ptr += alpha * (t_row_base * t_col_base);
660 ++t_col_base;
661 ++a_mat_ptr;
662 }
663 ++t_row_base;
664 }
665 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
666 ++t_row_base;
667 ++t_coords;
668 ++t_w; // move to another integration weight
669 }
670
671 OpBase::locMat *= vol;
672
673 MoFEMFunctionReturn(0);
674};
675
676template <int FIELD_DIM, typename OpBase>
677MoFEMErrorCode OpMassImpl<1, FIELD_DIM, GAUSS, OpBase>::integrateImpl(
678 EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data,
679 double vol) {
680 MoFEMFunctionBegin;
681
682 auto integrate = [&](auto &mat) {
683 MoFEMFunctionBeginHot;
684 // get integration weights
685 auto t_w = OpBase::getFTensor0IntegrationWeight();
686 // get base function gradient on rows
687 auto t_row_base = row_data.getFTensor0N();
688 // get coordinate at integration points
689 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
690
691 FTensor::Index<'i', FIELD_DIM> i;
692 auto get_t_vec = [&](const int rr) {
693 std::array<double *, FIELD_DIM> ptrs;
694 for (auto i = 0; i != FIELD_DIM; ++i)
695 ptrs[i] = &mat(rr + i, i);
696 return FTensor::Tensor1<FTensor::PackPtr<double *, FIELD_DIM>, FIELD_DIM>(
697 ptrs);
698 };
699
700 // loop over integration points
701 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
702 const double beta = betaCoeff(t_coords(0), t_coords(1), t_coords(2));
703 // take into account Jacobian
704 const double alpha = t_w * beta;
705 // loop over rows base functions
706 int rr = 0;
707 for (; rr != OpBase::nbRows / FIELD_DIM; rr++) {
708 // get column base functions gradient at gauss point gg
709 auto t_col_base = col_data.getFTensor0N(gg, 0);
710 // get mat vec
711 auto t_vec = get_t_vec(FIELD_DIM * rr);
712 // loop over columns
713 for (int cc = 0; cc != OpBase::nbCols / FIELD_DIM; cc++) {
714 // calculate element of local matrix
715 t_vec(i) += alpha * (t_row_base * t_col_base);
716 ++t_col_base;
717 ++t_vec;
718 }
719 ++t_row_base;
720 }
721 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
722 ++t_row_base;
723 ++t_coords;
724 ++t_w; // move to another integration weight
725 }
726 MoFEMFunctionReturnHot(0);
727 };
728
729 CHKERR integrate(OpBase::locMat);
730 OpBase::locMat *= vol;
731
732 MoFEMFunctionReturn(0);
733};
734
735template <int FIELD_DIM, typename OpBase>
736MoFEMErrorCode OpMassImpl<3, FIELD_DIM, GAUSS, OpBase>::iNtegrate(
737 EntitiesFieldData::EntData &row_data,
738 EntitiesFieldData::EntData &col_data) {
739 FTensor::Index<'i', FIELD_DIM> i;
740 MoFEMFunctionBegin;
741 size_t nb_base_functions = row_data.getN().size2() / 3;
742 // // get element volume
743 const double vol = OpBase::getMeasure();
744 // get integration weights
745 auto t_w = OpBase::getFTensor0IntegrationWeight();
746 // get base function gradient on rows
747 auto t_row_base = row_data.getFTensor1N<3>();
748 // get coordinate at integration points
749 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
750 // loop over integration points
751 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
752 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
753 // take into account Jacobian
754 const double alpha = t_w * beta;
755 // loop over rows base functions
756 auto a_mat_ptr = &*OpBase::locMat.data().begin();
757 int rr = 0;
758 for (; rr != OpBase::nbRows; rr++) {
759 // get column base functions gradient at gauss point gg
760 auto t_col_base = col_data.getFTensor1N<3>(gg, 0);
761 // loop over columns
762 for (int cc = 0; cc != OpBase::nbCols; cc++) {
763 // calculate element of local matrix
764 (*a_mat_ptr) += alpha * (t_row_base(i) * t_col_base(i));
765 ++t_col_base;
766 ++a_mat_ptr;
767 }
768 ++t_row_base;
769 }
770 for (; rr < nb_base_functions; ++rr)
771 ++t_row_base;
772 ++t_coords;
773 ++t_w; // move to another integration weight
774 }
775 MoFEMFunctionReturn(0);
776};
777
778template <typename OpBase>
779MoFEMErrorCode OpMassImpl<3, 4, GAUSS, OpBase>::iNtegrate(
780 EntitiesFieldData::EntData &row_data,
781 EntitiesFieldData::EntData &col_data) {
782 MoFEMFunctionBegin;
783 FTensor::Index<'i', 2> I;
784 FTensor::Index<'k', 3> k;
785 auto get_t_vec = [&](const int rr) {
786 return FTensor::Tensor1<FTensor::PackPtr<double *, 2>, 2>{
787 &OpBase::locMat(rr + 0, 0), &OpBase::locMat(rr + 1, 1)};
788 };
789 size_t nb_base_functions = row_data.getN().size2() / 3;
790 // // get element volume
791 const double vol = OpBase::getMeasure();
792 // get integration weights
793 auto t_w = OpBase::getFTensor0IntegrationWeight();
794 // get base function gradient on rows
795 auto t_row_base = row_data.getFTensor1N<3>();
796 // get coordinate at integration points
797 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
798 // loop over integration points
799 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
800 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
801 // take into account Jacobian
802 const double alpha = t_w * beta;
803 // loop over rows base functions
804 int rr = 0;
805 for (; rr != OpBase::nbRows / 2; rr++) {
806 // get column base functions gradient at gauss point gg
807 auto t_col_base = col_data.getFTensor1N<3>(gg, 0);
808 auto t_vec = get_t_vec(2 * rr);
809 // loop over columns
810 for (int cc = 0; cc != OpBase::nbCols / 2; cc++) {
811 // calculate element of local matrix
812 t_vec(I) += alpha * (t_row_base(k) * t_col_base(k));
813 ++t_col_base;
814 ++t_vec;
815 }
816 ++t_row_base;
817 }
818 for (; rr < nb_base_functions; ++rr)
819 ++t_row_base;
820 ++t_coords;
821 ++t_w; // move to another integration weight
822 }
823 MoFEMFunctionReturn(0);
824}
825
826template <typename OpBase>
827MoFEMErrorCode OpMassImpl<3, 9, GAUSS, OpBase>::iNtegrate(
828 EntitiesFieldData::EntData &row_data,
829 EntitiesFieldData::EntData &col_data) {
830 MoFEMFunctionBegin;
831 FTensor::Index<'i', 3> i;
832 FTensor::Index<'k', 3> k;
833 auto get_t_vec = [&](const int rr) {
834 return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>{
835 &OpBase::locMat(rr + 0, 0), &OpBase::locMat(rr + 1, 1),
836 &OpBase::locMat(rr + 2, 2)};
837 };
838 size_t nb_base_functions = row_data.getN().size2() / 3;
839 // // get element volume
840 const double vol = OpBase::getMeasure();
841 // get integration weights
842 auto t_w = OpBase::getFTensor0IntegrationWeight();
843 // get base function gradient on rows
844 auto t_row_base = row_data.getFTensor1N<3>();
845 // get coordinate at integration points
846 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
847 // loop over integration points
848 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
849 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
850 // take into account Jacobian
851 const double alpha = t_w * beta;
852 // loop over rows base functions
853 int rr = 0;
854 for (; rr != OpBase::nbRows / 3; rr++) {
855 // get column base functions gradient at gauss point gg
856 auto t_col_base = col_data.getFTensor1N<3>(gg, 0);
857 auto t_vec = get_t_vec(3 * rr);
858 // loop over columns
859 for (int cc = 0; cc != OpBase::nbCols / 3; cc++) {
860 // calculate element of local matrix
861 t_vec(i) += alpha * (t_row_base(k) * t_col_base(k));
862 ++t_col_base;
863 ++t_vec;
864 }
865 ++t_row_base;
866 }
867 for (; rr < nb_base_functions; ++rr)
868 ++t_row_base;
869 ++t_coords;
870 ++t_w; // move to another integration weight
871 }
872 MoFEMFunctionReturn(0);
873};
874
875template <typename OpBase>
876MoFEMErrorCode OpMassImpl<9, 9, GAUSS, OpBase>::iNtegrate(
877 EntitiesFieldData::EntData &row_data,
878 EntitiesFieldData::EntData &col_data) {
879 MoFEMFunctionBegin;
880 FTensor::Index<'i', 3> i;
881 FTensor::Index<'j', 3> j;
882 size_t nb_base_functions = row_data.getN().size2() / 9;
883 // // get element volume
884 const double vol = OpBase::getMeasure();
885 // get integration weights
886 auto t_w = OpBase::getFTensor0IntegrationWeight();
887 // get base function gradient on rows
888 auto t_row_base = row_data.getFTensor2N<3, 3>();
889 // get coordinate at integration points
890 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
891 // loop over integration points
892 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
893 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
894 // take into account Jacobian
895 const double alpha = t_w * beta;
896 // loop over rows base functions
897 int rr = 0;
898 for (; rr != OpBase::nbRows; rr++) {
899 // get column base functions gradient at gauss point gg
900 auto t_col_base = col_data.getFTensor2N<3, 3>(gg, 0);
901 // loop over columns
902 for (int cc = 0; cc != OpBase::nbCols; cc++) {
903 // calculate element of local matrix
904 OpBase::locMat(rr, cc) += alpha * (t_row_base(i, j) * t_col_base(i, j));
905 ++t_col_base;
906 }
907 ++t_row_base;
908 }
909 for (; rr < nb_base_functions; ++rr)
910 ++t_row_base;
911 ++t_coords;
912 ++t_w; // move to another integration weight
913 }
914 MoFEMFunctionReturn(0);
915}
916
917template <int FIELD_DIM, IntegrationType I, typename OpBase>
918MoFEMErrorCode OpMassCacheImpl<1, FIELD_DIM, I, OpBase>::iNtegrate(
919 EntitiesFieldData::EntData &row_data,
920 EntitiesFieldData::EntData &col_data) {
921 MoFEMFunctionBegin;
922
923 const auto vol = this->getMeasure();
924 const auto row_type = this->rowType;
925 const auto col_type = this->colType;
926 auto &loc_mat = this->locMat;
927
928 auto p = std::make_pair(row_type, col_type);
929
930 if (cacheLocMats[p]) {
931 if (cacheLocMats[p]->size1() != loc_mat.size1() &&
932 cacheLocMats[p]->size2() != loc_mat.size2()) {
933 cacheLocMats[p]->resize(loc_mat.size1(), loc_mat.size2());
934 CHKERR this->integrateImpl(row_data, col_data, 1);
935 *(cacheLocMats[p]) = loc_mat;
936 } else {
937 loc_mat = *(cacheLocMats[p]);
938 }
939 loc_mat *= scalarBeta * this->getMeasure();
940 } else {
941 CHKERR this->integrateImpl(row_data, col_data,
942 scalarBeta * this->getMeasure());
943 }
944 MoFEMFunctionReturn(0);
945}
946
947template <int SPACE_DIM, int S, typename OpBase>
948MoFEMErrorCode
949OpGradSymTensorGradImpl<1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase>::iNtegrate(
950 EntitiesFieldData::EntData &row_data,
951 EntitiesFieldData::EntData &col_data) {
952 MoFEMFunctionBegin;
953
954 const size_t nb_row_dofs = row_data.getIndices().size();
955 const size_t nb_col_dofs = col_data.getIndices().size();
956
957 if (nb_row_dofs && nb_col_dofs) {
958
959 const bool diag = (row_data.getFieldEntities()[0]->getLocalUniqueId() ==
960 col_data.getFieldEntities()[0]->getLocalUniqueId());
961
962 FTensor::Index<'i', SPACE_DIM> i;
963 FTensor::Index<'j', SPACE_DIM> j;
964 FTensor::Index<'k', SPACE_DIM> k;
965 FTensor::Index<'l', SPACE_DIM> l;
966
967 // get element volume
968 double vol = OpBase::getMeasure();
969
970 // get intergrayion weights
971 auto t_w = OpBase::getFTensor0IntegrationWeight();
972
973 // get derivatives of base functions on rows
974 auto t_row_diff_base = row_data.getFTensor1DiffN<SPACE_DIM>();
975
976 // Elastic stiffness tensor (4th rank tensor with minor and major
977 // symmetry)
978 using DL = DataLayoutTraits<DataLayout::GaussByCoeffs>;
979 auto get_D_at_pts =
980 MatrixSizeHelper<GetFTensor4DdgFromMatType<SPACE_DIM, SPACE_DIM, S, DL>,
981 DL>::get(*matD, OpBase::nbIntegrationPts);
982 auto t_D_at_pts = get_D_at_pts();
983
984 // get coordinate at integration points
985 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
986
987 // iterate over integration points
988 for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
989
990 // calculate scalar weight times element volume
991 double a = t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
992
993 // iterate over row base functions
994 int rr = 0;
995 for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
996
997 // get sub matrix for the row
998 auto t_m = OpBase::template getLocMat<SPACE_DIM>(SPACE_DIM * rr);
999
1000 FTensor::Christof<double, SPACE_DIM, SPACE_DIM> t_rowD;
1001 // I mix up the indices here so that it behaves like a
1002 // Dg. That way I don't have to have a separate wrapper
1003 // class Christof_Expr, which simplifies things.
1004 t_rowD(l, j, k) = t_D_at_pts(i, j, k, l) * (a * t_row_diff_base(i));
1005
1006 // get derivatives of base functions for columns
1007 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1008
1009 // iterate column base functions
1010 const auto nb_cols = (diag) ? rr : OpBase::nbCols / SPACE_DIM - 1;
1011 for (int cc = 0; cc <= nb_cols; ++cc) {
1012
1013 // integrate block local stiffens matrix
1014 t_m(i, j) += t_rowD(i, j, k) * t_col_diff_base(k);
1015
1016 // move to next column base function
1017 ++t_col_diff_base;
1018
1019 // move to next block of local stiffens matrix
1020 ++t_m;
1021 }
1022
1023 // move to next row base function
1024 ++t_row_diff_base;
1025 }
1026
1027 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1028 ++t_row_diff_base;
1029
1030 // move to next integration weight
1031 ++t_w;
1032 ++t_D_at_pts;
1033 ++t_coords;
1034 }
1035
1036 // Copy symmetry
1037 if (diag) {
1038 for (int rr = 0; rr != OpBase::nbRows / SPACE_DIM - 1; ++rr) {
1039 auto t_m_rr = getFTensor2FromArray<SPACE_DIM, SPACE_DIM, SPACE_DIM>(
1040 this->locMat, SPACE_DIM * rr, SPACE_DIM * (rr + 1));
1041 auto t_m_cc = getFTensor2FromArray<SPACE_DIM, SPACE_DIM>(
1042 this->locMat, SPACE_DIM * (rr + 1), SPACE_DIM * rr,
1043 SPACE_DIM * OpBase::nbCols);
1044 for (int cc = rr + 1; cc < OpBase::nbCols / SPACE_DIM; ++cc) {
1045 t_m_rr(i, j) = t_m_cc(j, i);
1046 ++t_m_rr;
1047 ++t_m_cc;
1048 }
1049 }
1050 }
1051 }
1052
1053 MoFEMFunctionReturn(0);
1054}
1055
1056template <int FIELD_DIM, int SPACE_DIM, int S, typename OpBase>
1057MoFEMErrorCode
1058OpGradGradSymTensorGradGradImpl<1, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase>::
1059 iNtegrate(EntitiesFieldData::EntData &row_data,
1060 EntitiesFieldData::EntData &col_data) {
1061 MoFEMFunctionBegin;
1062
1063 if (OpBase::nbRows && OpBase::nbCols) {
1064
1065 FTensor::Index<'i', SPACE_DIM> i;
1066 FTensor::Index<'j', SPACE_DIM> j;
1067 FTensor::Index<'k', SPACE_DIM> k;
1068 FTensor::Index<'l', SPACE_DIM> l;
1069 FTensor::Index<'I', FIELD_DIM> I;
1070 FTensor::Index<'J', FIELD_DIM> J;
1071
1072 auto &row_hessian = row_data.getN(BaseDerivatives::SecondDerivative);
1073 auto &col_hessian = col_data.getN(BaseDerivatives::SecondDerivative);
1074
1075#ifndef NDEBUG
1076 if (row_hessian.size1() != OpBase::nbIntegrationPts) {
1077 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1078 "Wrong number of integration pts (%zu != %d)", row_hessian.size1(),
1079 OpBase::nbIntegrationPts);
1080 }
1081 if (row_hessian.size2() !=
1082 OpBase::nbRowBaseFunctions * SPACE_DIM * SPACE_DIM) {
1083 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1084 "Wrong number of base functions (%zu != %d)",
1085 row_hessian.size2() / (SPACE_DIM * SPACE_DIM),
1086 OpBase::nbRowBaseFunctions);
1087 }
1088 if (row_hessian.size2() < OpBase::nbRows * SPACE_DIM * SPACE_DIM) {
1089 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1090 "Wrong number of base functions (%zu < %d)", row_hessian.size2(),
1091 OpBase::nbRows * SPACE_DIM * SPACE_DIM);
1092 }
1093 if (col_hessian.size1() != OpBase::nbIntegrationPts) {
1094 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1095 "Wrong number of integration pts (%zu != %d)", col_hessian.size1(),
1096 OpBase::nbIntegrationPts);
1097 }
1098 if (col_hessian.size2() < OpBase::nbCols * SPACE_DIM * SPACE_DIM) {
1099 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1100 "Wrong number of base functions (%zu < %d)", col_hessian.size2(),
1101 OpBase::nbCols * SPACE_DIM * SPACE_DIM);
1102 }
1103#endif
1104
1105 // get element volume
1106 double vol = OpBase::getMeasure();
1107
1108 // get intergrayion weights
1109 auto t_w = OpBase::getFTensor0IntegrationWeight();
1110
1111 auto t_row_diff2 = getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(
1112 &*row_hessian.data().begin());
1113
1114 // Elastic stiffness tensor (4th rank tensor with minor and major
1115 // symmetry)
1116 using DL = DataLayoutTraits<DataLayout::GaussByCoeffs>;
1117 auto get_D_at_pts =
1118 MatrixSizeHelper<GetFTensor4DdgFromMatType<SPACE_DIM, SPACE_DIM, S, DL>,
1119 DL>::get(*matD, OpBase::nbIntegrationPts);
1120 auto t_D_at_pts = get_D_at_pts();
1121 auto t_kd = FTensor::Kronecker_Delta<double>();
1122
1123 // iterate over integration points
1124 for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1125
1126 // calculate scalar weight times element volume
1127 double a = t_w * vol;
1128
1129 // iterate over row base functions
1130 int rr = 0;
1131 for (; rr != OpBase::nbRows; ++rr) {
1132
1133 auto t_mat = getFTensor2FromPtr<FIELD_DIM, FIELD_DIM, FIELD_DIM>(
1134 &OpBase::locMat(rr, 0));
1135
1136 FTensor::Tensor2_symmetric<double, SPACE_DIM> t_rowD;
1137 t_rowD(k, l) = t_D_at_pts(i, j, k, l) * (a * t_row_diff2(i, j));
1138
1139 // get derivatives of base functions for columns
1140 auto t_col_diff2 =
1141 getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(&col_hessian(gg, 0));
1142
1143 // iterate column base functions
1144 for (int cc = 0; cc != OpBase::nbCols; ++cc) {
1145
1146 // integrate block local stiffens matrix
1147 t_mat(I, J) += (t_rowD(i, j) * t_col_diff2(i, j)) * t_kd(I, J);
1148
1149 // move to next column base function
1150 ++t_col_diff2;
1151
1152 // move to next block of local stiffens matrix
1153 ++t_mat;
1154 }
1155
1156 // move to next row base function
1157 ++t_row_diff2;
1158 }
1159
1160 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1161 ++t_row_diff2;
1162
1163 // move to next integration weight
1164 ++t_w;
1165 ++t_D_at_pts;
1166 }
1167 }
1168
1169 MoFEMFunctionReturn(0);
1170}
1171
1172template <int FIELD_DIM, int SPACE_DIM, int S, typename OpBase>
1173MoFEMErrorCode
1174OpGradGradTensorGradGradImpl<1, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase>::
1175 iNtegrate(EntitiesFieldData::EntData &row_data,
1176 EntitiesFieldData::EntData &col_data) {
1177 MoFEMFunctionBegin;
1178
1179 if (OpBase::nbRows && OpBase::nbCols) {
1180
1181 FTensor::Index<'i', SPACE_DIM> i;
1182 FTensor::Index<'j', SPACE_DIM> j;
1183 FTensor::Index<'k', SPACE_DIM> k;
1184 FTensor::Index<'l', SPACE_DIM> l;
1185 FTensor::Index<'I', FIELD_DIM> I;
1186 FTensor::Index<'J', FIELD_DIM> J;
1187
1188 auto &row_hessian = row_data.getN(BaseDerivatives::SecondDerivative);
1189 auto &col_hessian = col_data.getN(BaseDerivatives::SecondDerivative);
1190
1191#ifndef NDEBUG
1192 if (row_hessian.size1() != OpBase::nbIntegrationPts) {
1193 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1194 "Wrong number of integration pts (%zu != %d)", row_hessian.size1(),
1195 OpBase::nbIntegrationPts);
1196 }
1197 if (row_hessian.size2() !=
1198 OpBase::nbRowBaseFunctions * SPACE_DIM * SPACE_DIM) {
1199 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1200 "Wrong number of base functions (%zu != %d)",
1201 row_hessian.size2() / (SPACE_DIM * SPACE_DIM),
1202 OpBase::nbRowBaseFunctions);
1203 }
1204 if (row_hessian.size2() < OpBase::nbRows * SPACE_DIM * SPACE_DIM) {
1205 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1206 "Wrong number of base functions (%zu < %d)", row_hessian.size2(),
1207 OpBase::nbRows * SPACE_DIM * SPACE_DIM);
1208 }
1209 if (col_hessian.size1() != OpBase::nbIntegrationPts) {
1210 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1211 "Wrong number of integration pts (%zu != %d)", col_hessian.size1(),
1212 OpBase::nbIntegrationPts);
1213 }
1214 if (col_hessian.size2() < OpBase::nbCols * SPACE_DIM * SPACE_DIM) {
1215 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1216 "Wrong number of base functions (%zu < %d)", col_hessian.size2(),
1217 OpBase::nbCols * SPACE_DIM * SPACE_DIM);
1218 }
1219#endif
1220
1221 // get element volume
1222 double vol = OpBase::getMeasure();
1223
1224 // get intergrayion weights
1225 auto t_w = OpBase::getFTensor0IntegrationWeight();
1226
1227 auto t_row_diff2 = getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(
1228 &*row_hessian.data().begin());
1229
1230 // Elastic stiffness tensor
1231 using DL = DataLayoutTraits<DataLayout::GaussByCoeffs>;
1232 auto get_D_at_pts = MatrixSizeHelper<
1233 GetFTensor4FromMatType<SPACE_DIM, SPACE_DIM, SPACE_DIM, SPACE_DIM, S,
1234 DL>,
1235 DL>::get(*matD, OpBase::nbIntegrationPts);
1236 auto t_D_at_pts = get_D_at_pts();
1237
1238 constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
1239
1240 // iterate over integration points
1241 for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1242
1243 // calculate scalar weight times element volume
1244 double a = t_w * vol;
1245
1246 // get sub matrix for the row
1247 auto m_ptr = &*OpBase::locMat.data().begin();
1248
1249 // iterate over row base functions
1250 int rr = 0;
1251 for (; rr != OpBase::nbRows; ++rr) {
1252
1253 FTensor::Tensor2_symmetric<double, SPACE_DIM> t_rowD;
1254 t_rowD(k, l) = t_D_at_pts(i, j, k, l) * (a * t_row_diff2(i, j));
1255
1256 // get derivatives of base functions for columns
1257 auto t_col_diff2 =
1258 getFTensor2SymmetricLowerFromPtr<SPACE_DIM>(&col_hessian(gg, 0));
1259
1260 auto t_mat = getFTensor2FromPtr<FIELD_DIM, FIELD_DIM, FIELD_DIM>(
1261 &OpBase::locMat(FIELD_DIM * rr, 0));
1262
1263 // iterate column base functions
1264 for (int cc = 0; cc != OpBase::nbCols; ++cc) {
1265
1266 // integrate block local stiffens matrix
1267 t_mat(I, J) += (t_rowD(i, j) * t_col_diff2(i, j)) * t_kd(I, J);
1268
1269 // move to next column base function
1270 ++t_col_diff2;
1271
1272 // move to next block of local stiffens matrix
1273 ++t_mat;
1274 }
1275
1276 // move to next row base function
1277 ++t_row_diff2;
1278 }
1279
1280 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1281 ++t_row_diff2;
1282
1283 // move to next integration weight
1284 ++t_w;
1285 ++t_D_at_pts;
1286 }
1287 }
1288
1289 MoFEMFunctionReturn(0);
1290}
1291
1292template <int SPACE_DIM, int S, typename OpBase>
1293MoFEMErrorCode
1294OpGradTensorGradImpl<1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase>::iNtegrate(
1295 EntitiesFieldData::EntData &row_data,
1296 EntitiesFieldData::EntData &col_data) {
1297 MoFEMFunctionBegin;
1298
1299 const size_t nb_row_dofs = row_data.getIndices().size();
1300 const size_t nb_col_dofs = col_data.getIndices().size();
1301
1302 if (nb_row_dofs && nb_col_dofs) {
1303
1304 FTensor::Index<'i', SPACE_DIM> i;
1305 FTensor::Index<'j', SPACE_DIM> j;
1306 FTensor::Index<'k', SPACE_DIM> k;
1307 FTensor::Index<'l', SPACE_DIM> l;
1308
1309 // get element volume
1310 double vol = OpBase::getMeasure();
1311
1312 // get intergrayion weights
1313 auto t_w = OpBase::getFTensor0IntegrationWeight();
1314
1315 // get derivatives of base functions on rows
1316 auto t_row_diff_base = row_data.getFTensor1DiffN<SPACE_DIM>();
1317
1318 // stiffness tensor (4th rank tensor)
1319 using DL = DataLayoutTraits<DataLayout::GaussByCoeffs>;
1320 auto get_D_at_pts = MatrixSizeHelper<
1321 GetFTensor4FromMatType<SPACE_DIM, SPACE_DIM, SPACE_DIM, SPACE_DIM, S,
1322 DL>,
1323 DL>::get(*matD, OpBase::nbIntegrationPts);
1324 auto t_D_at_pts = get_D_at_pts();
1325
1326 // get coordinate at integration points
1327 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1328
1329 // iterate over integration points
1330 for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1331
1332 // calculate scalar weight times element volume
1333 double a = t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1334
1335 // iterate over row base functions
1336 int rr = 0;
1337 for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1338
1339 // get sub matrix for the row
1340 auto t_m = OpBase::template getLocMat<SPACE_DIM>(SPACE_DIM * rr);
1341
1342 // calculate row
1343 FTensor::Tensor3<double, 3, 3, 3> t_row;
1344 t_row(i, k, l) = t_D_at_pts(i, j, k, l) * (a * t_row_diff_base(j));
1345
1346 // get derivatives of base functions for columns
1347 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1348
1349 // iterate column base functions
1350 for (int cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1351
1352 // integrate block local stiffens matrix
1353 t_m(i, k) += t_row(i, k, l) * t_col_diff_base(l);
1354
1355 // move to next column base function
1356 ++t_col_diff_base;
1357
1358 // move to next block of local stiffens matrix
1359 ++t_m;
1360 }
1361
1362 // move to next row base function
1363 ++t_row_diff_base;
1364 }
1365
1366 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1367 ++t_row_diff_base;
1368
1369 // move to next integration weight
1370 ++t_w;
1371 ++t_D_at_pts;
1372 ++t_coords;
1373 }
1374 }
1375
1376 MoFEMFunctionReturn(0);
1377}
1378
1379template <int SPACE_DIM, typename OpBase>
1380MoFEMErrorCode OpMixDivTimesScalarImpl<SPACE_DIM, GAUSS, OpBase>::iNtegrate(
1381 EntitiesFieldData::EntData &row_data,
1382 EntitiesFieldData::EntData &col_data) {
1383 MoFEMFunctionBegin;
1384
1385 auto t_w = this->getFTensor0IntegrationWeight();
1386
1387 size_t nb_base_functions = row_data.getN().size2() / 3;
1388 auto t_row_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
1389
1390 const double alpha_constant = alphaConstant();
1391
1392 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1393 const double alpha = alpha_constant * this->getMeasure() * t_w;
1394
1395 size_t rr = 0;
1396 for (; rr != OpBase::nbRows; ++rr) {
1397 const double t_row_div_base = t_row_diff_base(i, i);
1398 auto t_col_base = col_data.getFTensor0N(gg, 0);
1399 for (size_t cc = 0; cc != OpBase::nbCols; ++cc) {
1400 this->locMat(rr, cc) += alpha * t_row_div_base * t_col_base;
1401 ++t_col_base;
1402 }
1403 ++t_row_diff_base;
1404 }
1405 for (; rr < nb_base_functions; ++rr)
1406 ++t_row_diff_base;
1407
1408 ++t_w;
1409 }
1410
1411 MoFEMFunctionReturn(0);
1412}
1413
1414template <int SPACE_DIM, typename OpBase, CoordinateTypes CoordSys>
1415MoFEMErrorCode
1416OpMixDivTimesVecImpl<SPACE_DIM, GAUSS, OpBase, CoordSys>::iNtegrate(
1417 EntitiesFieldData::EntData &row_data,
1418 EntitiesFieldData::EntData &col_data) {
1419 MoFEMFunctionBegin;
1420
1421 auto t_w = this->getFTensor0IntegrationWeight();
1422 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1423
1424 size_t nb_base_functions = row_data.getN().size2() / 3;
1425 auto t_row_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
1426 auto t_row_base = row_data.getFTensor1N<3>();
1427 const double alpha_constant = alphaConstant() * this->getMeasure();
1428 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1429
1430 const double alpha =
1431 alpha_constant * t_w * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1432
1433 size_t rr = 0;
1434 for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1435 auto t_mat_diag = getFTensor1FromArrayDiag<SPACE_DIM, SPACE_DIM>(
1436 this->locMat, SPACE_DIM * rr);
1437 const double t_row_div_base = t_row_diff_base(i, i);
1438 auto t_col_base = col_data.getFTensor0N(gg, 0);
1439
1440 for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1441 t_mat_diag(i) += alpha * t_row_div_base * t_col_base;
1442 if constexpr (CoordSys == CYLINDRICAL) {
1443 t_mat_diag(i) += t_row_base(0) * (alpha / t_coords(0)) * t_col_base;
1444 }
1445 ++t_col_base;
1446 ++t_mat_diag;
1447 }
1448 ++t_row_base;
1449 ++t_row_diff_base;
1450 }
1451 for (; rr < nb_base_functions; ++rr) {
1452 ++t_row_diff_base;
1453 ++t_row_base;
1454 }
1455
1456 ++t_w;
1457 ++t_coords;
1458 }
1459
1460 MoFEMFunctionReturn(0);
1461}
1462
1463template <int SPACE_DIM, typename OpBase, CoordinateTypes COORDINATE_SYSTEM>
1464MoFEMErrorCode
1465OpMixScalarTimesDivImpl<SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM>::iNtegrate(
1466 EntitiesFieldData::EntData &row_data,
1467 EntitiesFieldData::EntData &col_data) {
1468 MoFEMFunctionBegin;
1469
1470#ifndef NDEBUG
1471 if (OpBase::locMat.size2() % SPACE_DIM)
1472 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1473 "Number of rows in matrix should be multiple of space dimensions");
1474#endif
1475
1476 // When we move to C++17 add if constexpr()
1477 if constexpr (COORDINATE_SYSTEM == POLAR || COORDINATE_SYSTEM == SPHERICAL)
1478 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1479 "%s coordiante not implemented",
1480 CoordinateTypesNames[COORDINATE_SYSTEM]);
1481
1482 auto t_w = this->getFTensor0IntegrationWeight();
1483 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1484 size_t nb_base_functions_row = row_data.getN().size2();
1485 auto t_row_base = row_data.getFTensor0N();
1486 const double vol = this->getMeasure();
1487 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1488
1489 const double alpha =
1490 alphaConstant(t_coords(0), t_coords(1), t_coords(2)) * t_w * vol;
1491
1492 size_t rr = 0;
1493 auto t_m = getFTensor1FromPtr<SPACE_DIM>(OpBase::locMat.data().data());
1494
1495 // When we move to C++17 add if constexpr()
1496 if constexpr (COORDINATE_SYSTEM == CARTESIAN) {
1497 for (; rr != OpBase::nbRows; ++rr) {
1498 const double r_val = alpha * t_row_base;
1499 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1500 for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1501 t_m(i) += r_val * t_col_diff_base(i);
1502 ++t_col_diff_base;
1503 ++t_m;
1504 }
1505 ++t_row_base;
1506 }
1507 }
1508
1509 // When we move to C++17 add if constexpr()
1510 if constexpr (COORDINATE_SYSTEM == CYLINDRICAL) {
1511 for (; rr != OpBase::nbRows; ++rr) {
1512 const double r_val = alpha * t_row_base;
1513 auto t_col_base = col_data.getFTensor0N(gg, 0);
1514 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1515 for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1516 t_m(i) += r_val * t_col_diff_base(i);
1517 t_m(0) += (r_val / t_coords(0)) * t_col_base;
1518 ++t_col_base;
1519 ++t_col_diff_base;
1520 ++t_m;
1521 }
1522 ++t_row_base;
1523 }
1524 }
1525
1526 for (; rr < nb_base_functions_row; ++rr)
1527 ++t_row_base;
1528
1529 ++t_w;
1530 ++t_coords;
1531 }
1532
1533 MoFEMFunctionReturn(0);
1534}
1535
1536template <int SPACE_DIM, typename OpBase>
1537MoFEMErrorCode
1538OpMixVectorTimesGradImpl<3, SPACE_DIM, SPACE_DIM, GAUSS, OpBase>::iNtegrate(
1539 EntitiesFieldData::EntData &row_data,
1540 EntitiesFieldData::EntData &col_data) {
1541 MoFEMFunctionBegin;
1542
1543 auto t_w = this->getFTensor0IntegrationWeight();
1544
1545 size_t nb_base_functions = row_data.getN().size2() / 3;
1546 auto t_row_base = row_data.getFTensor1N<3>();
1547 auto &mat = this->locMat;
1548 const double alpha_constant = alphaConstant();
1549 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1550
1551 const double alpha = alpha_constant * this->getMeasure() * t_w;
1552
1553 size_t rr = 0;
1554 for (; rr != OpBase::nbRows; ++rr) {
1555 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1556 for (size_t cc = 0; cc != OpBase::nbCols; ++cc) {
1557 mat(rr, cc) += alpha * t_row_base(i) * t_col_diff_base(i);
1558 ++t_col_diff_base;
1559 }
1560 ++t_row_base;
1561 }
1562 for (; rr < nb_base_functions; ++rr)
1563 ++t_row_base;
1564
1565 ++t_w;
1566 }
1567
1568 MoFEMFunctionReturn(0);
1569}
1570
1571template <int SPACE_DIM, typename OpBase>
1572MoFEMErrorCode
1573OpMixVectorTimesGradImpl<1, SPACE_DIM, SPACE_DIM, GAUSS, OpBase>::iNtegrate(
1574 EntitiesFieldData::EntData &row_data,
1575 EntitiesFieldData::EntData &col_data) {
1576 MoFEMFunctionBegin;
1577
1578 auto t_w = this->getFTensor0IntegrationWeight();
1579
1580 size_t nb_base_functions = row_data.getN().size2();
1581 auto t_row_base = row_data.getFTensor0N();
1582
1583 auto get_t_vec = [&](const int rr) {
1584 std::array<double *, SPACE_DIM> ptrs;
1585 for (auto i = 0; i != SPACE_DIM; ++i)
1586 ptrs[i] = &OpBase::locMat(rr + i, 0);
1587 return FTensor::Tensor1<FTensor::PackPtr<double *, 1>, SPACE_DIM>(ptrs);
1588 };
1589
1590 const double alpha_constant = alphaConstant();
1591 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1592
1593 const double alpha = alpha_constant * this->getMeasure() * t_w;
1594
1595 size_t rr = 0;
1596 for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1597 auto t_vec = get_t_vec(SPACE_DIM * rr);
1598 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1599 for (size_t cc = 0; cc != OpBase::nbCols; ++cc) {
1600 t_vec(i) += alpha * t_row_base * t_col_diff_base(i);
1601 ++t_col_diff_base;
1602 ++t_vec;
1603 }
1604 ++t_row_base;
1605 }
1606 for (; rr < nb_base_functions; ++rr)
1607 ++t_row_base;
1608
1609 ++t_w;
1610 }
1611
1612 MoFEMFunctionReturn(0);
1613}
1614
1615template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, typename OpBase>
1616MoFEMErrorCode
1617OpBaseTensorBaseImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase>::iNtegrate(
1618 EntitiesFieldData::EntData &row_data,
1619 EntitiesFieldData::EntData &col_data) {
1620 MoFEMFunctionBegin;
1621 constexpr bool tensor_basis = BASE_DIM == 4 || BASE_DIM == 9;
1622 constexpr int base_extent = tensor_basis ? (BASE_DIM == 4 ? 2 : 3) : BASE_DIM;
1623 constexpr int col_dim = tensor_basis ? 1 : SPACE_DIM;
1624 static_assert(base_extent >= SPACE_DIM,
1625 "Column basis storage is smaller than the physical tensor");
1626 using DL = DataLayoutTraits<DataLayout::GaussByCoeffs>;
1627 FTENSOR_INDEXES(SPACE_DIM, i, j);
1628 FTensor::Index<'a', FIELD_DIM> a;
1629 auto get_tensor = MatrixSizeHelper<
1630 GetFTensor3FromMatType<SPACE_DIM, SPACE_DIM, FIELD_DIM, S, DL>, DL>::get(
1631 *matTensor, OpBase::nbIntegrationPts);
1632 auto t_tensor = get_tensor();
1633 auto t_w = this->getFTensor0IntegrationWeight();
1634 auto t_row = row_data.getFTensor0N();
1635 const double alpha_constant = alphaConstant() * this->getMeasure();
1636
1637 for (int gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1638 const double alpha = alpha_constant * t_w;
1639 int rr = 0;
1640 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
1641 auto t_col = [&]() {
1642 if constexpr (tensor_basis)
1643 return col_data.getFTensor2N<base_extent, base_extent>(gg, 0);
1644 else
1645 return col_data.getFTensor1N<BASE_DIM>(gg, 0);
1646 }();
1647 auto t_m = OpBase::template getLocMat<FIELD_DIM, col_dim>(FIELD_DIM * rr);
1648 for (int cc = 0; cc != OpBase::nbCols / col_dim; ++cc) {
1649 if constexpr (tensor_basis) {
1650 const FTensor::Number<0> n0;
1651 t_m(a, n0) += (alpha * t_row) * (t_tensor(i, j, a) * t_col(i, j));
1652 } else {
1653 t_m(a, i) += (alpha * t_row) * (t_tensor(i, j, a) * t_col(j));
1654 }
1655 ++t_m;
1656 ++t_col;
1657 }
1658 ++t_row;
1659 }
1660 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1661 ++t_row;
1662 ++t_tensor;
1663 ++t_w;
1664 }
1665 MoFEMFunctionReturn(0);
1666}
1667
1668template <int SPACE_DIM, typename OpBase>
1669MoFEMErrorCode OpMixTensorTimesGradImpl<SPACE_DIM, GAUSS, OpBase>::iNtegrate(
1670 EntitiesFieldData::EntData &row_data,
1671 EntitiesFieldData::EntData &col_data) {
1672 MoFEMFunctionBegin;
1673
1674 auto t_w = this->getFTensor0IntegrationWeight();
1675 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1676
1677 size_t nb_base_functions = row_data.getN().size2() / 3;
1678 auto t_row_base = row_data.getFTensor1N<3>();
1679 const double alpha_constant = alphaConstant() * this->getMeasure();
1680 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1681
1682 const double alpha =
1683 alpha_constant * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * t_w;
1684
1685 size_t rr = 0;
1686 for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1687 auto t_mat_diag = getFTensor1FromArrayDiag<SPACE_DIM, SPACE_DIM>(
1688 this->locMat, SPACE_DIM * rr);
1689 auto t_col_diff_base = col_data.getFTensor1DiffN<SPACE_DIM>(gg, 0);
1690
1691 for (size_t cc = 0; cc != OpBase::nbCols / SPACE_DIM; ++cc) {
1692 t_mat_diag(i) += alpha * t_row_base(j) * t_col_diff_base(j);
1693 ++t_col_diff_base;
1694 ++t_mat_diag;
1695 }
1696
1697 ++t_row_base;
1698 }
1699 for (; rr < nb_base_functions; ++rr)
1700 ++t_row_base;
1701
1702 ++t_w;
1703 ++t_coords;
1704 }
1705
1706 MoFEMFunctionReturn(0);
1707}
1708
1709} // namespace MoFEM
1710
1711#endif //__BILINEAR_FORMS_INTEGRATORS_IMPL_HPP__
FTENSOR_INDEXES
#define FTENSOR_INDEXES(DIM,...)
Definition Templates.hpp:2672
a
constexpr double a
Definition approx_sphere.cpp:30
SPACE_DIM
constexpr int SPACE_DIM
Definition child_and_parent.cpp:17
FIELD_DIM
constexpr int FIELD_DIM
Definition child_and_parent.cpp:16
FTensor::Kronecker_Delta
Kronecker Delta class.
Definition Kronecker_Delta.hpp:15
FTensor::Tensor2_symmetric
Definition Tensor2_symmetric_value.hpp:14
CoordinateTypesNames
static const char *const CoordinateTypesNames[]
Coordinate system names.
Definition definitions.h:121
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition definitions.h:463
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition definitions.h:359
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32
CoordinateTypes
CoordinateTypes
Coodinate system.
Definition definitions.h:127
CYLINDRICAL
@ CYLINDRICAL
Definition definitions.h:130
POLAR
@ POLAR
Definition definitions.h:129
CARTESIAN
@ CARTESIAN
Definition definitions.h:128
SPHERICAL
@ SPHERICAL
Definition definitions.h:131
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition definitions.h:432
CHKERR
#define CHKERR
Inline error check.
Definition definitions.h:551
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition definitions.h:456
BASE_DIM
constexpr int BASE_DIM
Definition dg_projection.cpp:35
t_kd
constexpr auto t_kd
Definition free_surface.cpp:161
MoFEM::IntegrationType
IntegrationType
Form integrator integration types.
Definition FormsIntegrators.hpp:237
MoFEM::ScalarFun
boost::function< double(const double, const double, const double)> ScalarFun
Scalar function type.
Definition FormsIntegrators.hpp:249
MoFEM::scalar_fun_one
double scalar_fun_one(const double, const double, const double)
Default scalar function returning 1.
Definition FormsIntegrators.hpp:260
MoFEM::GAUSS
@ GAUSS
Gaussian quadrature integration.
Definition FormsIntegrators.hpp:238
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition LogManager.hpp:316
i
FTensor::Index< 'i', SPACE_DIM > i
Definition hcurl_divergence_operator_2d.cpp:27
J
FTensor::Index< 'J', DIM1 > J
Definition level_set.cpp:30
l
FTensor::Index< 'l', 3 > l
Definition matrix_function.cpp:23
j
FTensor::Index< 'j', 3 > j
Definition matrix_function.cpp:21
k
FTensor::Index< 'k', 3 > k
Definition matrix_function.cpp:22
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition Exceptions.hpp:56
MoFEM
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
MoFEM::GetFTensor4FromMatType
decltype(GetFTensor4FromMatImpl< Tensor_Dim0, Tensor_Dim1, Tensor_Dim2, Tensor_Dim3, S, DL, M >::get(std::declval< M & >(), 0, 0)) GetFTensor4FromMatType
Definition Templates.hpp:1355
MoFEM::CacheMatsTypeType
std::map< std::pair< int, int >, boost::shared_ptr< MatrixDouble > > CacheMatsTypeType
Definition BiLinearFormsIntegratorsImpl.hpp:19
MoFEM::ConstantFun
boost::function< double()> ConstantFun
Constant function type.
Definition FormsIntegrators.hpp:280
MoFEM::GetFTensor3FromMatType
decltype(GetFTensor3FromMatImpl< Tensor_Dim0, Tensor_Dim1, Tensor_Dim2, S, DL, M >::get(std::declval< M & >(), 0, 0)) GetFTensor3FromMatType
Definition Templates.hpp:1008
I
constexpr IntegrationType I
Definition operators_tests.cpp:31
OpBase
OpBaseImpl< PETSC, EdgeEleOp > OpBase
Definition radiation.cpp:29
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition EntitiesFieldData.hpp:150
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:297
MoFEM::EntitiesFieldData::EntData::getFTensor0N
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
Definition EntitiesFieldData.hpp:1806
MoFEM::EntitiesFieldData::EntData::getFieldEntities
const VectorFieldEntities & getFieldEntities() const
Get field entities (const version)
Definition EntitiesFieldData.hpp:1591
MoFEM::EntitiesFieldData::EntData::getFTensor2N
auto getFTensor2N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
Definition EntitiesFieldData.hpp:1139
MoFEM::EntitiesFieldData::EntData::getFTensor1DiffN
auto getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
Definition EntitiesFieldData.hpp:802
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:1600
MoFEM::EntitiesFieldData::EntData::getFTensor1N
auto getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
Definition EntitiesFieldData.hpp:948
MoFEM::EntitiesFieldData::EntData::getIndices
const VectorInt & getIndices() const
Get global indices of degrees of freedom on entity.
Definition EntitiesFieldData.hpp:1523
MoFEM::OpBaseImpl::rowSide
int rowSide
row side number
Definition FormsIntegrators.hpp:374
MoFEM::OpBaseImpl::colSide
int colSide
column side number
Definition FormsIntegrators.hpp:375
MoFEM::OpBaseImpl::nbRows
int nbRows
number of dofs on rows
Definition FormsIntegrators.hpp:369
MoFEM::OpBaseImpl::colType
EntityType colType
column type
Definition FormsIntegrators.hpp:377
MoFEM::OpBaseImpl::nbIntegrationPts
int nbIntegrationPts
number of integration points
Definition FormsIntegrators.hpp:371
MoFEM::OpBaseImpl::locMat
MatrixDouble locMat
local entity block matrix
Definition FormsIntegrators.hpp:382
MoFEM::OpBaseImpl::nbCols
int nbCols
number if dof on column
Definition FormsIntegrators.hpp:370
MoFEM::OpBaseImpl::nbRowBaseFunctions
int nbRowBaseFunctions
number or row base functions
Definition FormsIntegrators.hpp:372
MoFEM::OpBaseImpl::rowType
EntityType rowType
row type
Definition FormsIntegrators.hpp:376
MoFEM::OpBaseTensorBaseImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase >::matTensor
boost::shared_ptr< MatrixDouble > matTensor
Definition BiLinearFormsIntegratorsImpl.hpp:454
MoFEM::OpBaseTensorBaseImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase >::alphaConstant
ConstantFun alphaConstant
Definition BiLinearFormsIntegratorsImpl.hpp:455
MoFEM::OpBaseTensorBaseImpl< BASE_DIM, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase >::OpBaseTensorBaseImpl
OpBaseTensorBaseImpl(const std::string row_field_name, const std::string col_field_name, boost::shared_ptr< MatrixDouble > mat_tensor, ConstantFun alpha_fun=[]() constexpr { return 1.;}, const bool assemble_transpose=false, const bool only_transpose=false)
Definition BiLinearFormsIntegratorsImpl.hpp:441
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::OpGradGradImpl< 1, 1, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:27
MoFEM::OpGradGradImpl< 1, 1, SPACE_DIM, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition BiLinearFormsIntegratorsImpl.hpp:39
MoFEM::OpGradGradImpl< 1, FIELD_DIM, SPACE_DIM, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition BiLinearFormsIntegratorsImpl.hpp:58
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, FIELD_DIM, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:46
MoFEM::OpGradGradSymTensorGradGradImpl< 1, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase >::matD
boost::shared_ptr< MatrixDouble > matD
Definition BiLinearFormsIntegratorsImpl.hpp:269
MoFEM::OpGradGradSymTensorGradGradImpl< 1, FIELD_DIM, 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:258
MoFEM::OpGradGradTensorGradGradImpl< 1, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase >::matD
boost::shared_ptr< MatrixDouble > matD
Definition BiLinearFormsIntegratorsImpl.hpp:292
MoFEM::OpGradGradTensorGradGradImpl< 1, FIELD_DIM, SPACE_DIM, S, GAUSS, OpBase >::OpGradGradTensorGradGradImpl
OpGradGradTensorGradGradImpl(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:281
MoFEM::OpGradSymTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition BiLinearFormsIntegratorsImpl.hpp:224
MoFEM::OpGradSymTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:210
MoFEM::OpGradSymTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::matD
boost::shared_ptr< MatrixDouble > matD
Definition BiLinearFormsIntegratorsImpl.hpp:223
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:211
MoFEM::OpGradTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::matD
boost::shared_ptr< MatrixDouble > matD
Definition BiLinearFormsIntegratorsImpl.hpp:245
MoFEM::OpGradTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition BiLinearFormsIntegratorsImpl.hpp:246
MoFEM::OpGradTensorGradImpl< 1, SPACE_DIM, SPACE_DIM, S, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:236
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:237
MoFEM::OpMassCacheImpl< 1, FIELD_DIM, I, OpBase >::scalarBeta
double scalarBeta
Definition BiLinearFormsIntegratorsImpl.hpp:199
MoFEM::OpMassCacheImpl< 1, FIELD_DIM, I, OpBase >::cacheLocMats
CacheMatsTypeType cacheLocMats
Definition BiLinearFormsIntegratorsImpl.hpp:200
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:187
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::OpMassImpl< 1, 1, GAUSS, OpBase >::iNtegrate
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
Integrate grad-grad operator.
Definition BiLinearFormsIntegratorsImpl.hpp:87
MoFEM::OpMassImpl< 1, FIELD_DIM, GAUSS, OpBase >::iNtegrate
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &data)
Class dedicated to integrate operator.
Definition BiLinearFormsIntegratorsImpl.hpp:107
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:103
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:132
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:149
MoFEM::OpMassImpl< 3, FIELD_DIM, GAUSS, OpBase >::betaCoeff
ScalarFun betaCoeff
Definition BiLinearFormsIntegratorsImpl.hpp:125
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:115
MoFEM::OpMassImpl< 9, 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:166
MoFEM::OpMixDivTimesScalarImpl< SPACE_DIM, GAUSS, OpBase >::alphaConstant
ConstantFun alphaConstant
Definition BiLinearFormsIntegratorsImpl.hpp:315
MoFEM::OpMixDivTimesScalarImpl< SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:314
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:302
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:328
MoFEM::OpMixDivTimesVecImpl< SPACE_DIM, GAUSS, OpBase, CoordSys >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:349
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:338
MoFEM::OpMixScalarTimesDivImpl< SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM >::alphaConstant
ScalarFun alphaConstant
Definition BiLinearFormsIntegratorsImpl.hpp:378
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:365
MoFEM::OpMixScalarTimesDivImpl< SPACE_DIM, GAUSS, OpBase, COORDINATE_SYSTEM >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:377
MoFEM::OpMixTensorTimesGradImpl< SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:488
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:466
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:476
MoFEM::OpMixTensorTimesGradImpl< SPACE_DIM, GAUSS, OpBase >::j
FTensor::Index< 'j', SPACE_DIM > j
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:489
MoFEM::OpMixVectorTimesGradImpl< 1, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::alphaConstant
ConstantFun alphaConstant
Definition BiLinearFormsIntegratorsImpl.hpp:424
MoFEM::OpMixVectorTimesGradImpl< 1, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:423
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:411
MoFEM::OpMixVectorTimesGradImpl< 3, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::i
FTensor::Index< 'i', SPACE_DIM > i
summit Index
Definition BiLinearFormsIntegratorsImpl.hpp:402
MoFEM::OpMixVectorTimesGradImpl< 3, SPACE_DIM, SPACE_DIM, GAUSS, OpBase >::alphaConstant
ConstantFun alphaConstant
Definition BiLinearFormsIntegratorsImpl.hpp:403
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:390
Generated by Doxygen 1.9.8 and hosted at University of Glasgow