v0.13.1
LinearFormsIntegrators.hpp
Go to the documentation of this file.
1/** \file LinearFormsIntegrators.hpp
2 * \brief Linear forms integrators
3 * \ingroup mofem_form
4
5*/
6
7
8
9#ifndef __LINEAR_FORMS_INTEGRATORS_HPP__
10#define __LINEAR_FORMS_INTEGRATORS_HPP__
11
12namespace MoFEM {
13
14template <int BASE_DIM, int FIELD_DIM, IntegrationType I, typename OpBase>
16
17/**
18 * @brief Integrate source
19 *
20 * @tparam OpBase
21 */
22template <typename OpBase>
23struct OpSourceImpl<1, 1, GAUSS, OpBase> : public OpBase {
24
25 /**
26 * @brief Construct a new Op Source Impl object
27 *
28 * @param field_name
29 * @param source_fun
30 * @param ents_ptr
31 */
32 OpSourceImpl(const std::string field_name, ScalarFun source_fun,
33 boost::shared_ptr<Range> ents_ptr = nullptr)
34 : OpBase(field_name, field_name, OpBase::OPROW), sourceFun(source_fun),
35 entsPtr(ents_ptr) {}
36
37protected:
38 boost::shared_ptr<Range> entsPtr;
41};
42
43template <int FIELD_DIM, typename OpBase>
44struct OpSourceImpl<1, FIELD_DIM, GAUSS, OpBase> : public OpBase {
45 OpSourceImpl(const std::string field_name, VectorFun<FIELD_DIM> source_fun,
46 boost::shared_ptr<Range> ents_ptr = nullptr)
47 : OpBase(field_name, field_name, OpBase::OPROW), sourceFun(source_fun),
48 entsPtr(ents_ptr) {}
49
50protected:
51 boost::shared_ptr<Range> entsPtr;
54};
55
56template <int BASE_DIM, typename OpBase>
58 OpSourceImpl(const std::string field_name, VectorFun<BASE_DIM> source_fun,
59 boost::shared_ptr<Range> ents_ptr = nullptr)
60 : OpBase(field_name, field_name, OpBase::OPROW), sourceFun(source_fun),
61 entsPtr(ents_ptr) {}
62
63protected:
64 boost::shared_ptr<Range> entsPtr;
67};
68
69template <int BASE_DIM, int S, IntegrationType I, typename OpBase>
71
72template <int S, typename OpBase>
74
76 boost::shared_ptr<VectorDouble> vec,
77 ScalarFun beta_coeff,
78 boost::shared_ptr<Range> ents_ptr = nullptr)
79 : OpBase(field_name, field_name, OpBase::OPROW), sourceVec(vec),
80 betaCoeff(beta_coeff), entsPtr(ents_ptr) {}
81
82protected:
84 boost::shared_ptr<Range> entsPtr;
85 boost::shared_ptr<VectorDouble> sourceVec;
87};
88
89template <int BASE_DIM, int FIELD_DIM, int S, IntegrationType I,
90 typename OpBase>
92
93template <int FIELD_DIM, int S, typename OpBase>
95
97 boost::shared_ptr<MatrixDouble> vec,
98 ScalarFun beta_coeff,
99 boost::shared_ptr<Range> ents_ptr = nullptr)
100 : OpBase(field_name, field_name, OpBase::OPROW), sourceVec(vec),
101 betaCoeff(beta_coeff), entsPtr(ents_ptr) {}
102
103protected:
105 boost::shared_ptr<Range> entsPtr;
106 boost::shared_ptr<MatrixDouble> sourceVec;
109};
110
111template <int BASE_DIM, int S, typename OpBase>
113 : public OpBase {
114
116 boost::shared_ptr<MatrixDouble> vec,
117 ScalarFun beta_coeff,
118 boost::shared_ptr<Range> ents_ptr = nullptr)
119 : OpBase(field_name, field_name, OpBase::OPROW), sourceVec(vec),
120 betaCoeff(beta_coeff), entsPtr(ents_ptr) {}
121
122protected:
124 boost::shared_ptr<Range> entsPtr;
125 boost::shared_ptr<MatrixDouble> sourceVec;
128};
129
130template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
131 typename OpBase>
133
134template <int SPACE_DIM, int S, typename OpBase>
136 : public OpBase {
137
139
141 boost::shared_ptr<MatrixDouble> mat_vals,
142 ScalarFun beta_coeff,
143 boost::shared_ptr<Range> ents_ptr = nullptr)
144 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
145 betaCoeff(beta_coeff) {}
146
147protected:
148 boost::shared_ptr<MatrixDouble> matVals;
149 boost::shared_ptr<Range> entsPtr;
152};
153
154template <int SPACE_DIM, int S, typename OpBase>
156 : public OpBase {
157
160
162 const std::string field_name, boost::shared_ptr<MatrixDouble> mat_vals,
163 ScalarFun beta_coeff = [](double, double, double) { return 1; },
164 boost::shared_ptr<Range> ents_ptr = nullptr)
165 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
166 betaCoeff(beta_coeff) {}
167
168protected:
169 boost::shared_ptr<MatrixDouble> matVals;
170 boost::shared_ptr<Range> entsPtr;
173};
174
175template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
176 typename OpBase>
178
179template <int SPACE_DIM, int S, typename OpBase>
181 : public OpBase {
182
184 const std::string field_name, boost::shared_ptr<MatrixDouble> &mat_vals,
185 ScalarFun beta_coeff = [](double, double, double) { return 1; })
186 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
187 betaCoeff(beta_coeff) {}
188
189protected:
190 boost::shared_ptr<MatrixDouble> matVals;
195};
196
197template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
198 typename OpBase>
200
201template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
203 : public OpBase {
205 boost::shared_ptr<MatrixDouble> mat_vals, ScalarFun beta,
206 boost::shared_ptr<Range> ents_ptr = nullptr)
207 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
208 betaConst(beta), entsPtr(ents_ptr) {}
209
210protected:
212 boost::shared_ptr<Range> entsPtr;
213 boost::shared_ptr<MatrixDouble> matVals;
217};
218
219template <int SPACE_DIM, typename OpBase>
222 boost::shared_ptr<VectorDouble> vec_vals, ScalarFun beta,
223 boost::shared_ptr<Range> ents_ptr = nullptr)
224 : OpBase(field_name, field_name, OpBase::OPROW), vecVals(vec_vals),
225 betaConst(beta), entsPtr(ents_ptr) {}
226
227protected:
229 boost::shared_ptr<Range> entsPtr;
230 boost::shared_ptr<VectorDouble> vecVals;
233};
234
235template <int FIELD_DIM, typename OpBase>
237 : public OpBase {
238
240 boost::shared_ptr<VectorDouble> vec, ScalarFun beta,
241 boost::shared_ptr<Range> ents_ptr = nullptr)
242 : OpBase(field_name, field_name, OpBase::OPROW), sourceVec(vec),
243 betaCoeff(beta), entsPtr(ents_ptr) {}
244
245protected:
247 boost::shared_ptr<Range> entsPtr;
248 boost::shared_ptr<VectorDouble> sourceVec;
250};
251
252template <int SPACE_DIM, IntegrationType I, typename OpBase>
254
255template <int SPACE_DIM, typename OpBase>
258 boost::shared_ptr<MatrixDouble> &mat_vals)
259 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals) {}
260
261protected:
262 boost::shared_ptr<MatrixDouble> matVals;
265};
266
267template <int SPACE_DIM, IntegrationType I, typename OpBase>
269
270template <int SPACE_DIM, typename OpBase>
273 boost::shared_ptr<MatrixDouble> &mat_vals)
274 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals) {}
275
276protected:
277 boost::shared_ptr<MatrixDouble> matVals;
281};
282
283/**
284 * @brief Multiply vactor times normal on the face times scalar function
285 *
286 * This operator typically will be used to evaluate natural boundary conditions
287 * for mixed formulation.
288 *
289 * @tparam BASE_DIM
290 * @tparam SPACE_DIM
291 * @tparam OpBase
292 */
293template <int SPACE_DIM, IntegrationType I, typename OpBase>
295
296template <typename OpBase>
299 ScalarFun source_fun,
300 boost::shared_ptr<Range> ents_ptr = nullptr)
301 : OpBase(field_name, field_name, OpBase::OPROW), sourceFun(source_fun),
302 entsPtr(ents_ptr) {}
303
304protected:
305 boost::shared_ptr<Range> entsPtr;
309};
310
311template <typename OpBase>
314 ScalarFun source_fun,
315 boost::shared_ptr<Range> ents_ptr = nullptr)
316 : OpBase(field_name, field_name, OpBase::OPROW), sourceFun(source_fun),
317 entsPtr(ents_ptr) {}
318
319protected:
320 boost::shared_ptr<Range> entsPtr;
324};
325
326template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
327 typename OpBase>
329
330template <int SPACE_DIM, typename OpBase>
333 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
334 boost::shared_ptr<MatrixDouble> y_grad_ptr,
335 ConstantFun source_fun = []() { return 1; })
336 : OpBase(field_name, field_name, OpBase::OPROW), uPtr(u_ptr),
337 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
338
339protected:
340 boost::shared_ptr<MatrixDouble> uPtr;
341 boost::shared_ptr<MatrixDouble> yGradPtr;
344};
345
346template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
348 : public OpBase {
350 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
351 boost::shared_ptr<MatrixDouble> y_grad_ptr,
352 ConstantFun source_fun = []() { return 1; })
353 : OpBase(field_name, field_name, OpBase::OPROW), uPtr(u_ptr),
354 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
355
356protected:
357 boost::shared_ptr<MatrixDouble> uPtr;
358 boost::shared_ptr<MatrixDouble> yGradPtr;
361};
362
363/**
364 * @brief Linear integrator form
365 * @ingroup mofem_forms
366 *
367 * @tparam EleOp
368 * @tparam A
369 * @tparam I
370 */
371template <typename EleOp>
372template <AssemblyType A>
373template <IntegrationType I>
374struct FormsIntegrators<EleOp>::Assembly<A>::LinearForm {
375
376 /**
377 * @brief Integrate \f$(v,f(\mathbf{x}))_\Omega\f$, f is a scalar
378 * @ingroup mofem_forms
379 *
380 * @note \f$f(\mathbf{x})\$ is scalar function of coordinates
381 *
382 * @tparam BASE_DIM
383 * @tparam FIELD_DIM
384 */
385 template <int BASE_DIM, int FIELD_DIM>
386 struct OpSource : public OpSourceImpl<BASE_DIM, FIELD_DIM, I, OpBase> {
387 using OpSourceImpl<BASE_DIM, FIELD_DIM, I, OpBase>::OpSourceImpl;
388 };
389
390 /**
391 * @brief Vector field integrator \f$(v_i,f)_\Omega\f$, f is a vector
392 * @ingroup mofem_forms
393 *
394 * @tparam BASE_DIM
395 */
396 template <int BASE_DIM, int S = 1>
398 : public OpBaseTimesScalarFieldImpl<BASE_DIM, S, I, OpBase> {
399 using OpBaseTimesScalarFieldImpl<BASE_DIM, S, I,
400 OpBase>::OpBaseTimesScalarFieldImpl;
401 };
402
403 /**
404 * @brief Vector field integrator \f$(v,f_i)_\Omega\f$, f is a vector
405 * @ingroup mofem_forms
406 *
407 * @tparam BASE_DIM
408 * @tparam FIELD_DIM
409 * @tparam 0
410 */
411 template <int BASE_DIM, int FIELD_DIM, int S>
412 struct OpBaseTimesVector
413 : public OpBaseTimesVectorImpl<BASE_DIM, FIELD_DIM, S, I, OpBase> {
414 using OpBaseTimesVectorImpl<BASE_DIM, FIELD_DIM, S, I,
415 OpBase>::OpBaseTimesVectorImpl;
416 };
417
418 //! [Source operator]
419
420 //! [Grad times tensor]
421
422 /**
423 * @brief Integrate \f$(v_{,i},f_{ij})_\Omega\f$, f is a vector
424 * @ingroup mofem_forms
425 *
426 * @note \f$f_{ij}\$ is tensor at integration points
427 *
428 * @tparam BASE_DIM
429 * @tparam FIELD_DIM
430 * @tparam SPACE_DIM
431 */
432 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
433 struct OpGradTimesTensor
434 : public OpGradTimesTensorImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, I,
435 OpBase> {
437 OpBase>::OpGradTimesTensorImpl;
438 };
439
440 /**
441 * @brief Integrate \f$(v_{,i},f_{ij})_\Omega\f$, f is symmetric tensor
442 * @ingroup mofem_forms
443 *
444 * @note \f$f_{ij}\$ is tensor at integration points
445 *
446 * @tparam BASE_DIM
447 * @tparam FIELD_DIM
448 * @tparam SPACE_DIM
449 */
450 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
451 struct OpGradTimesSymTensor
452 : public OpGradTimesSymTensorImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, S, I,
453 OpBase> {
455 OpBase>::OpGradTimesSymTensorImpl;
456 };
457
458 /**
459 * @brief Integrate \f$(\lambda_{ij,j},u_{i})_\Omega\f$
460 *
461 * @tparam BASE_DIM
462 * @tparam FIELD_DIM
463 * @tparam SPACE_DIM
464 */
465 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
466 struct OpMixDivTimesU
467 : public OpMixDivTimesUImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, I, OpBase> {
469 OpBase>::OpMixDivTimesUImpl;
470 };
471
472 /**
473 * @brief Integrate \f$(\lambda_{ij},u_{i,j})_\Omega\f$
474 *
475 * @tparam SPACE_DIM
476 */
477 template <int SPACE_DIM>
478 struct OpMixTensorTimesGradU
479 : public OpMixTensorTimesGradUImpl<SPACE_DIM, I, OpBase> {
481 OpBase>::OpMixTensorTimesGradUImpl;
482 };
483
484 /**
485 * @brief Integrate \f$(u_{i},\lambda_{ij,j})_\Omega\f$
486 *
487 * @tparam SPACE_DIM
488 */
489 template <int SPACE_DIM>
490 struct OpMixVecTimesDivLambda
491 : public OpMixVecTimesDivLambdaImpl<SPACE_DIM, I, OpBase> {
492 using OpMixVecTimesDivLambdaImpl<SPACE_DIM, I,
493 OpBase>::OpMixVecTimesDivLambdaImpl;
494 };
495
496 /**
497 * @brief Multiply vactor times normal on the face times scalar function
498 *
499 * This operator typically will be used to evaluate natural boundary
500 * conditions for mixed formulation.
501 *
502 * @tparam BASE_DIM
503 * @tparam SPACE_DIM
504 * @tparam OpBase
505 */
506 template <int SPACE_DIM>
507 struct OpNormalMixVecTimesScalar
508 : public OpNormalMixVecTimesScalarImpl<SPACE_DIM, I, OpBase> {
509 using OpNormalMixVecTimesScalarImpl<SPACE_DIM, I,
510 OpBase>::OpNormalMixVecTimesScalarImpl;
511 };
512
513 /**
514 * @brief Convective term
515 *
516 * \f[
517 * (v, u_i \mathbf{y}_{,i})
518 * \f]
519 * where \f$\mathbf{y}\$ can be scalar or vector.
520 *
521 * @tparam BASE_DIM
522 * @tparam FIELD_DIM
523 * @tparam SPACE_DIM
524 */
525 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
526 struct OpConvectiveTermRhs
527 : public OpConvectiveTermRhsImpl<BASE_DIM, FIELD_DIM, SPACE_DIM, I,
528 OpBase> {
530 OpBase>::OpConvectiveTermRhsImpl;
531 };
532};
533
534template <typename OpBase>
536 EntitiesFieldData::EntData &row_data) {
538
539 if (entsPtr) {
540 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
542 }
543
544 // get element volume
545 const double vol = OpBase::getMeasure();
546 // get integration weights
547 auto t_w = OpBase::getFTensor0IntegrationWeight();
548 // get base function gradient on rows
549 auto t_row_base = row_data.getFTensor0N();
550 // get coordinate at integration points
551 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
552 // loop over integration points
553 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
554 // take into account Jacobian
555 const double alpha =
556 t_w * vol * sourceFun(t_coords(0), t_coords(1), t_coords(2));
557 // loop over rows base functions
558 int rr = 0;
559 for (; rr != OpBase::nbRows; ++rr) {
560 OpBase::locF[rr] += alpha * t_row_base;
561 ++t_row_base;
562 }
563 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
564 ++t_row_base;
565 ++t_coords;
566 ++t_w; // move to another integration weight
567 }
569}
570
571template <int FIELD_DIM, typename OpBase>
573 EntitiesFieldData::EntData &row_data) {
576 if (entsPtr) {
577 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
579 }
580
581 // get element volume
582 const double vol = OpBase::getMeasure();
583 // get integration weights
584 auto t_w = OpBase::getFTensor0IntegrationWeight();
585 // get base function gradient on rows
586 auto t_row_base = row_data.getFTensor0N();
587 // get coordinate at integration points
588 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
589 // loop over integration points
590 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
591 // source file
592 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
593 // take into account Jacobian
594 const double alpha = t_w * vol;
595 // loop over rows base functions
596 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(OpBase::locF);
597 int rr = 0;
598 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
599 t_nf(i) += alpha * t_row_base * t_source(i);
600 ++t_row_base;
601 ++t_nf;
602 }
603 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
604 ++t_row_base;
605 ++t_coords;
606 ++t_w; // move to another integration weight
607 }
609}
610
611template <int BASE_DIM, typename OpBase>
613 EntitiesFieldData::EntData &row_data) {
616 if (entsPtr) {
617 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
619 }
620
621 const size_t nb_base_functions = row_data.getN().size2() / BASE_DIM;
622 // get element volume
623 const double vol = OpBase::getMeasure();
624 // get integration weights
625 auto t_w = OpBase::getFTensor0IntegrationWeight();
626 // get base function gradient on rows
627 auto t_row_base = row_data.getFTensor1N<BASE_DIM>();
628 // get coordinate at integration points
629 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
630 // loop over integration points
631 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
632 // source file
633 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
634 // take into account Jacobian
635 const double alpha = t_w * vol;
636 // loop over rows base functions
637 int rr = 0;
638 for (; rr != OpBase::nbRows; ++rr) {
639 OpBase::locF[rr] += alpha * t_row_base(i) * t_source(i);
640 ++t_row_base;
641 }
642 for (; rr < nb_base_functions; ++rr)
643 ++t_row_base;
644 ++t_coords;
645 ++t_w; // move to another integration weight
646 }
648}
649
650template <int S, typename OpBase>
652 EntitiesFieldData::EntData &row_data) {
654 if (entsPtr) {
655 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
657 }
658 // get element volume
659 const double vol = OpBase::getMeasure();
660 // get integration weights
661 auto t_w = OpBase::getFTensor0IntegrationWeight();
662 // get base function gradient on rows
663 auto t_row_base = row_data.getFTensor0N();
664 // get vector values
665 auto t_vec = getFTensor0FromVec<S>(*sourceVec);
666
667#ifndef NDEBUG
668 if (sourceVec->size() != OpBase::nbIntegrationPts) {
669 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
670 "Wrong number of integration points %d != %d",
671 OpBase::nbIntegrationPts, sourceVec->size());
672 }
673#endif
674
675 // get coords
676 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
677 // loop over integration points
678 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
679 // take into account Jacobian
680 const double alpha =
681 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
682 // loop over rows base functions
683 int rr = 0;
684 for (; rr != OpBase::nbRows; ++rr) {
685 OpBase::locF[rr] += alpha * t_row_base * t_vec;
686 ++t_row_base;
687 }
688 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
689 ++t_row_base;
690 ++t_w; // move to another integration weight
691 ++t_vec;
692 ++t_coords;
693 }
695}
696
697template <int FIELD_DIM, int S, typename OpBase>
699 EntitiesFieldData::EntData &row_data) {
701 if (entsPtr) {
702 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
704 }
705 // get element volume
706 const double vol = OpBase::getMeasure();
707 // get integration weights
708 auto t_w = OpBase::getFTensor0IntegrationWeight();
709 // get base function gradient on rows
710 auto t_row_base = row_data.getFTensor0N();
711 // get coords
712 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
713 // get vector values
714 auto t_vec = getFTensor1FromMat<FIELD_DIM, S>(*sourceVec);
715 // loop over integration points
716 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
717 // take into account Jacobian
718 const double alpha =
719 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
720 // get loc vector tensor
721 auto t_nf = OpBase::template getNf<FIELD_DIM>();
722 // loop over rows base functions
723 int rr = 0;
724 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
725 t_nf(i) += alpha * t_row_base * t_vec(i);
726 ++t_row_base;
727 ++t_nf;
728 }
729 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
730 ++t_row_base;
731 ++t_w; // move to another integration weight
732 ++t_vec;
733 ++t_coords;
734 }
736}
737
738template <int BASE_DIM, int S, typename OpBase>
741 EntitiesFieldData::EntData &row_data) {
743 if (entsPtr) {
744 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
746 }
747 const size_t nb_base_functions = row_data.getN().size2() / BASE_DIM;
748 // get element volume
749 const double vol = OpBase::getMeasure();
750 // get integration weights
751 auto t_w = OpBase::getFTensor0IntegrationWeight();
752 // get base function gradient on rows
753 auto t_row_base = row_data.getFTensor1N<BASE_DIM>();
754 // get coords
755 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
756 // get vector values
757 auto t_vec = getFTensor1FromMat<BASE_DIM, S>(*sourceVec);
758 // loop over integration points
759 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
760 // take into account Jacobian
761 const double alpha =
762 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
763 // loop over rows base functions
764 int rr = 0;
765 for (; rr != OpBase::nbRows; ++rr) {
766 OpBase::locF[rr] += alpha * t_row_base(i) * t_vec(i);
767 ++t_row_base;
768 }
769 for (; rr < nb_base_functions; ++rr)
770 ++t_row_base;
771 ++t_w; // move to another integration weight
772 ++t_vec;
773 ++t_coords;
774 }
776}
777
778template <int SPACE_DIM, int S, typename OpBase>
781 EntitiesFieldData::EntData &row_data) {
783 if (entsPtr) {
784 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
786 }
787 // get element volume
788 const double vol = OpBase::getMeasure();
789 // get integration weights
790 auto t_w = OpBase::getFTensor0IntegrationWeight();
791 // get base function gradient on rows
792 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
793 // get filed gradient values
794 auto t_val_grad = getFTensor1FromMat<SPACE_DIM, S>(*(matVals));
795 // get coordinate at integration points
796 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
797 // loop over integration points
798 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
799 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
800 // take into account Jacobian
801 const double alpha = t_w * beta;
802 // loop over rows base functions
803 int rr = 0;
804 for (; rr != OpBase::nbRows; rr++) {
805 // calculate element of local matrix
806 OpBase::locF[rr] += alpha * (t_row_grad(i) * t_val_grad(i));
807 ++t_row_grad; // move to another element of gradient of base
808 // function on row
809 }
810 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
811 ++t_row_grad;
812
813 ++t_coords;
814 ++t_val_grad;
815 ++t_w; // move to another integration weight
816 }
818}
819
820template <int SPACE_DIM, int S, typename OpBase>
823 EntitiesFieldData::EntData &row_data) {
825 if (entsPtr) {
826 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
828 }
829 // get element volume
830 const double vol = OpBase::getMeasure();
831 // get integration weights
832 auto t_w = OpBase::getFTensor0IntegrationWeight();
833 // get base function gradient on rows
834 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
835 // get filed gradient values
836 auto t_val_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
837 // get coordinate at integration points
838 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
839 // loop over integration points
840 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
841 // take into account Jacobian
842 const double alpha =
843 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
844 // get rhs vector
845 auto t_nf = OpBase::template getNf<SPACE_DIM>();
846 // loop over rows base functions
847 int rr = 0;
848 for (; rr != OpBase::nbRows / SPACE_DIM; rr++) {
849 // calculate element of local matrix
850 t_nf(i) += alpha * (t_row_grad(j) * t_val_grad(i, j));
851 ++t_row_grad; // move to another element of gradient of base
852 // function on row
853 ++t_nf;
854 }
855 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
856 ++t_row_grad;
857
858 ++t_coords;
859 ++t_val_grad;
860 ++t_w; // move to another integration weight
861 }
863}
864
865template <int SPACE_DIM, int S, typename OpBase>
868 EntitiesFieldData::EntData &row_data) {
870 // get element volume
871 const double vol = OpBase::getMeasure();
872 // get coords
873 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
874 // get integration weights
875 auto t_w = OpBase::getFTensor0IntegrationWeight();
876 // get base function gradient on rows
877 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
878 // get filed gradient values
879 auto t_val_mat = getFTensor2SymmetricFromMat<SPACE_DIM, S>(*(matVals));
880 // loop over integration points
881 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
882 // take into account Jacobian
883 const double alpha =
884 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
885 // get rhs vector
886 auto t_nf = OpBase::template getNf<SPACE_DIM>();
887 // loop over rows base functions
888 int rr = 0;
889 for (; rr != OpBase::nbRows / SPACE_DIM; rr++) {
890 // calculate element of local matrix
891 t_nf(j) += alpha * (t_row_grad(i) * t_val_mat(i, j));
892 ++t_row_grad; // move to another element of gradient of base
893 // function on row
894 ++t_nf;
895 }
896 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
897 ++t_row_grad;
898 ++t_val_mat;
899 ++t_coords;
900 ++t_w; // move to another integration weight
901 }
903}
904
905template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
908 EntitiesFieldData::EntData &row_data) {
910
911 if (entsPtr) {
912 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
914 }
915
916 const size_t nb_base_functions = row_data.getN().size2() / 3;
917 auto t_w = this->getFTensor0IntegrationWeight();
918 // get coordinate at integration points
919 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
920 auto t_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
921 auto t_u = getFTensor1FromMat<FIELD_DIM>(*(matVals));
922
923 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
924
925 const double alpha = this->getMeasure() * t_w *
926 betaConst(t_coords(0), t_coords(1), t_coords(2));
927 auto t_nf = OpBase::template getNf<FIELD_DIM>();
928
929 size_t bb = 0;
930 for (; bb != this->nbRows / FIELD_DIM; ++bb) {
931 const double t_div_base = t_diff_base(j, j);
932 t_nf(i) += alpha * t_div_base * t_u(i);
933 ++t_nf;
934 ++t_diff_base;
935 }
936 for (; bb < nb_base_functions; ++bb)
937 ++t_diff_base;
938
939 ++t_u;
940 ++t_w;
941 ++t_coords;
942 }
943
945}
946
947template <int SPACE_DIM, typename OpBase>
949 EntitiesFieldData::EntData &row_data) {
951
952 if (entsPtr) {
953 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
955 }
956
957 const size_t nb_base_functions = row_data.getN().size2() / 3;
958 auto t_w = this->getFTensor0IntegrationWeight();
959 // get coordinate at integration points
960 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
961 auto t_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
962 auto t_u = getFTensor0FromVec(*(vecVals));
963
964 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
965
966 const double alpha = this->getMeasure() * t_w *
967 betaConst(t_coords(0), t_coords(1), t_coords(2));
968 ;
969
970 size_t bb = 0;
971 for (; bb != this->nbRows; ++bb) {
972 const double t_div_base = t_diff_base(j, j);
973 OpBase::locF[bb] += alpha * t_div_base * t_u;
974 ++t_diff_base;
975 }
976 for (; bb < nb_base_functions; ++bb)
977 ++t_diff_base;
978
979 ++t_u;
980 ++t_w;
981 ++t_coords;
982 }
983
985}
986
987/**
988 * @brief div U times vector
989 *
990 * \f[
991 * \delta u_j = \phi^m\delta\overline{u}^m_j\\
992 * \delta u_{j,i} = \phi^m_{,i}\delta\overline{u}^m_j\\
993 * \textrm{tr}[\delta u_{j,i}] = \delta u_{j,i}\delta_{ji}\\
994 * (\textrm{tr}[\delta u_{j,i}], v) =\\
995 * (\delta u_{j,i} \delta_{ij}, v) =\\
996 * (\delta u_{j,i}, \delta_{ij} v) =\\
997 * (\phi^m_{,i}\delta\overline{u}^m_j, \delta_{ij} v) \\
998 * f_i^m=(\phi^m_{,i}, v)
999 * \f]
1000 *
1001 * @tparam FIELD_DIM
1002 * @tparam SPACE_DIM
1003 * @tparam OpBase
1004 * @param row_data
1005 * @return MoFEMErrorCode
1006 */
1007template <int FIELD_DIM, typename OpBase>
1010 EntitiesFieldData::EntData &row_data) {
1013
1014 if (entsPtr) {
1015 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
1017 }
1018
1019 // get element volume
1020 const double vol = OpBase::getMeasure();
1021 // get integration weights
1022 auto t_w = OpBase::getFTensor0IntegrationWeight();
1023 // get coordinate at integration points
1024 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1025 // get base function gradient on rows
1026 auto t_row_grad = row_data.getFTensor1DiffN<FIELD_DIM>();
1027 // get vector values
1028 auto t_vec = getFTensor0FromVec(*sourceVec);
1029 // loop over integration points
1030 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1031 // take into account Jacobian
1032 const double alpha =
1033 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1034 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(OpBase::locF);
1035 // loop over rows base functions
1036 int rr = 0;
1037 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
1038 t_nf(i) += alpha * t_row_grad(i) * t_vec;
1039 ++t_row_grad;
1040 ++t_nf;
1041 }
1042 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1043 ++t_row_grad;
1044 ++t_w; // move to another integration weight
1045 ++t_vec;
1046 ++t_coords;
1047 }
1049}
1050
1051template <int SPACE_DIM, typename OpBase>
1053 EntitiesFieldData::EntData &row_data) {
1055
1056 const size_t nb_base_functions = row_data.getN().size2() / 3;
1057 auto t_w = this->getFTensor0IntegrationWeight();
1058 auto t_base = row_data.getFTensor1N<3>();
1059 auto t_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
1060
1061 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1062
1063 const double alpha = this->getMeasure() * t_w;
1064 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1065
1066 size_t bb = 0;
1067 for (; bb != this->nbRows / SPACE_DIM; ++bb) {
1068 t_nf(i) += alpha * t_base(j) * t_grad(i, j);
1069 ++t_nf;
1070 ++t_base;
1071 }
1072 for (; bb < nb_base_functions; ++bb)
1073 ++t_base;
1074
1075 ++t_grad;
1076 ++t_w;
1077 }
1078
1080}
1081
1082template <int SPACE_DIM, typename OpBase>
1084 EntitiesFieldData::EntData &row_data) {
1086
1087 const size_t nb_base_functions = row_data.getN().size2();
1088 auto t_w = this->getFTensor0IntegrationWeight();
1089 auto t_base = row_data.getFTensor0N();
1090 auto t_div = getFTensor1FromMat<SPACE_DIM>(*(matVals));
1091
1092 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1093 const double alpha = this->getMeasure() * t_w;
1094 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1095
1096 size_t bb = 0;
1097 for (; bb != this->nbRows / SPACE_DIM; ++bb) {
1098 t_nf(i) += alpha * t_base * t_div(i);
1099 ++t_nf;
1100 ++t_base;
1101 }
1102 for (; bb < nb_base_functions; ++bb)
1103 ++t_base;
1104
1105 ++t_div;
1106 ++t_w;
1107 }
1108
1110}
1111
1112template <typename OpBase>
1114 EntitiesFieldData::EntData &row_data) {
1116 if (entsPtr) {
1117 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
1119 }
1120 const size_t nb_base_functions = row_data.getN().size2() / 3;
1121 // get element volume
1122 const double vol = OpBase::getMeasure();
1123 // get integration weights
1124 auto t_w = OpBase::getFTensor0IntegrationWeight();
1125 // get base function gradient on rows
1126 auto t_row_base = row_data.getFTensor1N<3>();
1127 // get coordinate at integration points
1128 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1129 // get normal
1130 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1131 // loop over integration points
1132 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1133 // take into account Jacobian
1134 const double alpha = t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2));
1135 // loop over rows base functions
1136 int rr = 0;
1137 for (; rr != OpBase::nbRows; ++rr) {
1138 OpBase::locF[rr] += alpha * t_row_base(i) * t_normal(i);
1139 ++t_row_base;
1140 }
1141 for (; rr < nb_base_functions; ++rr)
1142 ++t_row_base;
1143 ++t_coords;
1144 ++t_w; // move to another integration weight
1145 ++t_normal;
1146 }
1148}
1149
1150template <typename OpBase>
1152 EntitiesFieldData::EntData &row_data) {
1154 if (entsPtr) {
1155 if (entsPtr->find(OpBase::getFEEntityHandle()) == entsPtr->end())
1157 }
1158 const size_t nb_base_functions = row_data.getN().size2() / 3;
1159 FTensor::Tensor1<double, 3> t_z{0., 0., 1.};
1160 // get element volume
1161 const double vol = OpBase::getMeasure();
1162 // get integration weights
1163 auto t_w = OpBase::getFTensor0IntegrationWeight();
1164 // get base function gradient on rows
1165 auto t_row_base = row_data.getFTensor1N<3>();
1166 // get coordinate at integration points
1167 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1168 // get normal
1169 auto t_tangent = OpBase::getFTensor1TangentAtGaussPts();
1170 // loop over integration points
1171 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1172 // take into account Jacobian
1173 const double alpha = t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2));
1177 t_normal(i) = FTensor::levi_civita(i, j, k) * t_tangent(j) * t_z(k);
1178 int rr = 0;
1179 for (; rr != OpBase::nbRows; ++rr) {
1180 OpBase::locF[rr] += alpha * t_row_base(i) * t_normal(i);
1181 ++t_row_base;
1182 }
1183 for (; rr < nb_base_functions; ++rr)
1184 ++t_row_base;
1185 ++t_coords;
1186 ++t_tangent;
1187 ++t_w; // move to another integration weight
1188 }
1190}
1191
1192template <int SPACE_DIM, typename OpBase>
1195 EntitiesFieldData::EntData &row_data) {
1197
1198 const size_t nb_base_functions = row_data.getN().size2();
1199 auto t_w = this->getFTensor0IntegrationWeight();
1200 auto t_base = row_data.getFTensor0N();
1201
1202 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1203 auto t_grad_y = getFTensor1FromMat<SPACE_DIM>(*yGradPtr);
1204
1206 const double alpha_constant = alphaConstant();
1207 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1208
1209 // get element volume
1210 const double vol = OpBase::getMeasure();
1211 const double c = (t_grad_y(i) * t_u(i)) * (t_w * vol * alpha_constant);
1212
1213 // get element volume
1214 int rr = 0;
1215 for (; rr != OpBase::nbRows; ++rr) {
1216 OpBase::locF[rr] += c * t_base;
1217 ++t_base;
1218 }
1219 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1220 ++t_base;
1221
1222 ++t_w; // move to another integration weight
1223 ++t_u;
1224 ++t_grad_y;
1225 }
1226
1228}
1229
1230template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
1233 EntitiesFieldData::EntData &row_data) {
1235
1236 const size_t nb_base_functions = row_data.getN().size2();
1237 auto t_w = this->getFTensor0IntegrationWeight();
1238 auto t_base = row_data.getFTensor0N();
1239
1240 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1241 auto t_grad_y = getFTensor2FromMat<FIELD_DIM, SPACE_DIM>(*yGradPtr);
1242
1245 const double alpha_constant = alphaConstant();
1246 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1247
1248 // get element volume
1249 const double vol = OpBase::getMeasure();
1250
1252 t_c(J) = (t_grad_y(J, i) * t_u(i)) * (t_w * vol * alpha_constant);
1253
1254 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(OpBase::locF);
1255 int rr = 0;
1256 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
1257 t_nf(J) += t_base * t_c(J);
1258 ++t_base;
1259 ++t_nf;
1260 }
1261 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1262 ++t_base;
1263 ++t_w; // move to another integration weight
1264
1265 ++t_u;
1266 ++t_grad_y;
1267 }
1268
1270}
1271
1272} // namespace MoFEM
1273
1274#endif // __LINEAR_FORMS_INTEGRATORS_HPP__
static Index< 'J', 3 > J
static Index< 'I', 3 > I
constexpr int SPACE_DIM
constexpr int FIELD_DIM
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
FormsIntegrators< DomainEleOp >::Assembly< PETSC >::LinearForm< GAUSS >::OpBaseTimesScalarField< 1, 1 > OpBaseTimesScalarField
IntegrationType
Form integrator integration types.
boost::function< double(const double, const double, const double)> ScalarFun
Scalar function type.
boost::function< FTensor::Tensor1< double, DIM >(const double, const double, const double)> VectorFun
Vector function type.
constexpr int BASE_DIM
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: MoFEM.hpp:24
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
boost::function< double()> ConstantFun
Constant function type.
constexpr auto field_name
OpBaseImpl< PETSC, EdgeEleOp > OpBase
Definition: radiation.cpp:29
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN(FieldApproximationBase base)
Get derivatives of base functions for Hdiv space.
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
VectorDouble locF
local entity vector
int nbRows
number of dofs on rows
int nbIntegrationPts
number of integration points
int nbRowBaseFunctions
number or row base functions
OpBaseTimesScalarFieldImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
OpConvectiveTermRhsImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > y_grad_ptr, ConstantFun source_fun=[]() { return 1;})
OpConvectiveTermRhsImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > y_grad_ptr, ConstantFun source_fun=[]() { return 1;})
OpGradTimesSymTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > &mat_vals, ScalarFun beta_coeff=[](double, double, double) { return 1;})
OpGradTimesTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff, boost::shared_ptr< Range > ents_ptr=nullptr)
OpGradTimesTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff=[](double, double, double) { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec_vals, ScalarFun beta, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixTensorTimesGradUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > &mat_vals)
OpMixVecTimesDivLambdaImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > &mat_vals)
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Multiply vactor times normal on the face times scalar function.
OpSourceImpl(const std::string field_name, ScalarFun source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Construct a new Op Source Impl object.
OpSourceImpl(const std::string field_name, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
OpSourceImpl(const std::string field_name, VectorFun< BASE_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)