v0.15.0
Loading...
Searching...
No Matches
LinearFormsIntegratorsImpl.hpp
Go to the documentation of this file.
1/** \file LinearFormsIntegratorsImpl.hpp
2 * \brief Linear forms integrators (implementation)
3 * \ingroup mofem_form
4
5*/
6
7#ifndef __LINEAR_FORMS_INTEGRATORS_IMPL_HPP__
8#define __LINEAR_FORMS_INTEGRATORS_IMPL_HPP__
9
10namespace MoFEM {
11
13 template <typename OpBase> struct S { S() = delete; };
15};
16
18 template <typename OpBase> struct S { S() = delete; };
20};
21
22template <int BASE_DIM, int FIELD_DIM, IntegrationType I, typename OpBase>
24
25/**
26 * @brief Integrate source
27 *
28 * @tparam OpBase
29 */
30template <typename OpBase>
32 : public OpBase {
33
34 /**
35 * @brief Construct a new Op Source Impl object
36 *
37 * @param field_name
38 * @param source_fun
39 * @param ents_ptr
40 */
41 OpSourceImpl(const std::string field_name, ScalarFun source_fun,
42 boost::shared_ptr<Range> ents_ptr = nullptr)
43 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
44 sourceFun(source_fun) {}
45
46 /**
47 * @brief Construct a new Op Source Impl object
48 *
49 * @param field_name
50 * @param time_fun
51 * @param source_fun
52 * @param ents_ptr
53 */
54 OpSourceImpl(const std::string field_name, TimeFun time_fun,
55 ScalarFun source_fun,
56 boost::shared_ptr<Range> ents_ptr = nullptr)
57 : OpBase(field_name, field_name, OpBase::OPROW, time_fun, ents_ptr),
58 sourceFun(source_fun) {}
59
60protected:
63};
64
65template <int FIELD_DIM, typename OpBase>
67 SourceFunctionSpecialization::S<OpBase>> : public OpBase {
68
69 /**
70 * @brief Construct a new Op Source Impl object
71 *
72 * @param field_name
73 * @param time_fun
74 * @param source_fun
75 * @param ents_ptr
76 */
77 OpSourceImpl(const std::string field_name, TimeFun time_fun,
78 VectorFun<FIELD_DIM> source_fun,
79 boost::shared_ptr<Range> ents_ptr = nullptr)
80 : OpBase(field_name, field_name, OpBase::OPROW, time_fun, ents_ptr) {}
81
82 /**
83 * @brief Construct a new Op Source Impl object
84 *
85 * @param field_name
86 * @param source_fun
87 * @param ents_ptr
88 */
89 OpSourceImpl(const std::string field_name, VectorFun<FIELD_DIM> source_fun,
90 boost::shared_ptr<Range> ents_ptr = nullptr)
91 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
92 sourceFun(source_fun) {}
93
94protected:
97};
98
99template <int FIELD_DIM, typename OpBase>
101 SourceFunctionSpecialization::S<OpBase>> : public OpBase {
102
103 OpSourceImpl(const std::string field_name, TimeFun time_fun,
104 VectorFun<FIELD_DIM> source_fun,
105 boost::shared_ptr<Range> ents_ptr = nullptr)
106 : OpBase(field_name, field_name, OpBase::OPROW, time_fun, ents_ptr),
107 sourceFun(source_fun) {}
108
109 OpSourceImpl(const std::string field_name, VectorFun<FIELD_DIM> source_fun,
110 boost::shared_ptr<Range> ents_ptr = nullptr)
111 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
112 sourceFun(source_fun) {}
113
114protected:
117};
118
119template <int BASE_DIM, int S, IntegrationType I, typename OpBase>
121
122template <int S, typename OpBase>
123struct OpBaseTimesScalarImpl<1, S, GAUSS, OpBase> : public OpBase {
124
126 const std::string field_name, boost::shared_ptr<VectorDouble> vec,
127 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; },
128 boost::shared_ptr<Range> ents_ptr = nullptr)
129 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr), sourceVec(vec),
130 betaCoeff(beta_coeff) {}
131
132protected:
134 boost::shared_ptr<VectorDouble> sourceVec;
136};
137
138template <int BASE_DIM, int FIELD_DIM, int S, IntegrationType I,
139 typename OpBase>
141
142template <int FIELD_DIM, int S, typename OpBase>
144
146 const std::string field_name, boost::shared_ptr<MatrixDouble> vec,
147 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; },
148 boost::shared_ptr<Range> ents_ptr = nullptr)
149 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr), sourceVec(vec),
150 betaCoeff(beta_coeff) {}
151
152protected:
154 boost::shared_ptr<MatrixDouble> sourceVec;
157};
158
159template <int FIELD_DIM, int S, typename OpBase>
161
163 const std::string field_name, boost::shared_ptr<MatrixDouble> vec,
164 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; },
165 boost::shared_ptr<Range> ents_ptr = nullptr)
166 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr), sourceVec(vec),
167 betaCoeff(beta_coeff) {}
168
169protected:
171 boost::shared_ptr<MatrixDouble> sourceVec;
174};
175
176template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
177 typename OpBase>
179
180template <int SPACE_DIM, int S, typename OpBase>
182 : public OpBase {
183
184 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
185
187 const std::string field_name, boost::shared_ptr<MatrixDouble> mat_vals,
188 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; },
189 boost::shared_ptr<Range> ents_ptr = nullptr)
190 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
191 matVals(mat_vals), betaCoeff(beta_coeff) {}
192
193protected:
194 boost::shared_ptr<MatrixDouble> matVals;
197};
198
199template <int SPACE_DIM, int S, typename OpBase>
201 : public OpBase {
202
203 FTensor::Index<'i', SPACE_DIM> i; ///< summit Index
204 FTensor::Index<'j', SPACE_DIM> j; ///< summit Index
205
207 const std::string field_name, boost::shared_ptr<MatrixDouble> mat_vals,
208 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; },
209 boost::shared_ptr<Range> ents_ptr = nullptr)
210 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
211 matVals(mat_vals), betaCoeff(beta_coeff) {}
212
213protected:
214 boost::shared_ptr<MatrixDouble> matVals;
217};
218
219template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S, IntegrationType I,
220 typename OpBase>
222
223template <int SPACE_DIM, int S, typename OpBase>
225 : public OpBase {
226
228 const std::string field_name, boost::shared_ptr<MatrixDouble> mat_vals,
229 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; })
230 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
231 betaCoeff(beta_coeff) {}
232
233protected:
234 boost::shared_ptr<MatrixDouble> matVals;
239};
240
241template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
242 typename OpBase, CoordinateTypes CoordSys>
244
245template <int FIELD_DIM, int SPACE_DIM, typename OpBase,
246 CoordinateTypes CoordSys>
248 : public OpBase {
250 const std::string field_name, boost::shared_ptr<MatrixDouble> mat_vals,
251 ScalarFun beta = [](double, double, double) { return 1; },
252 boost::shared_ptr<Range> ents_ptr = nullptr)
253 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
254 matVals(mat_vals), betaConst(beta) {}
255
256protected:
258 boost::shared_ptr<MatrixDouble> matVals;
262};
263
264template <int SPACE_DIM, typename OpBase, CoordinateTypes CoordSys>
265struct OpMixDivTimesUImpl<3, 1, SPACE_DIM, GAUSS, OpBase, CoordSys>
266 : public OpBase {
268 const std::string field_name, boost::shared_ptr<VectorDouble> vec_vals,
269 ScalarFun beta = [](double, double, double) constexpr { return 1; },
270 boost::shared_ptr<Range> ents_ptr = nullptr)
271 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
272 vecVals(vec_vals), betaConst(beta) {}
273
274protected:
276 boost::shared_ptr<VectorDouble> vecVals;
279};
280
281template <int FIELD_DIM, typename OpBase, CoordinateTypes CoordSys>
283 : public OpBase {
284
286 const std::string field_name, boost::shared_ptr<VectorDouble> vec,
287 ScalarFun beta = [](double, double, double) constexpr { return 1; },
288 boost::shared_ptr<Range> ents_ptr = nullptr)
289 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr), sourceVec(vec),
290 betaCoeff(beta) {}
291
292protected:
294 boost::shared_ptr<VectorDouble> sourceVec;
296};
297
298/**
299 * @brief Vector field time divergence of tensor
300 *
301 * \f[
302 * (v_i,\lambda_{ij,j})_\Omega
303 * \f]
304 *
305 * @tparam SPACE_DIM
306 * @tparam I
307 * @tparam OpBase
308 */
309template <int SPACE_DIM, IntegrationType I, typename OpBase>
311
312template <int SPACE_DIM, typename OpBase>
315 boost::shared_ptr<MatrixDouble> mat_vals)
316 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals) {}
317
319 boost::shared_ptr<MatrixDouble> mat_vals,
320 ScalarFun beta_fun)
321 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
322 betaCoeff(beta_fun) {}
323
324protected:
325 boost::shared_ptr<MatrixDouble> matVals;
326 ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
327
330};
331
332/**
333 * @brief Tensor field time gradient of vector field
334 *
335 * \f[
336 * (u_i,\lambda_{ij,j})_\Omega
337 * \f]
338 *
339 * @tparam SPACE_DIM
340 * @tparam I
341 * @tparam OpBase
342 */
343template <int SPACE_DIM, IntegrationType I, typename OpBase>
345
346template <int SPACE_DIM, typename OpBase>
349 boost::shared_ptr<MatrixDouble> mat_vals)
350 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals) {}
351
353 boost::shared_ptr<MatrixDouble> mat_vals,
354 ScalarFun beta_fun)
355 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
356 betaCoeff(beta_fun) {}
357
358protected:
359 boost::shared_ptr<MatrixDouble> matVals;
360 ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
361
365};
366
367/**
368 * @brief Multiply vector times normal on the face times scalar function
369 *
370 * This operator typically will be used to evaluate natural boundary conditions
371 * for mixed formulation.
372 *
373 * @tparam BASE_DIM
374 * @tparam SPACE_DIM
375 * @tparam OpBase
376 */
377template <int SPACE_DIM, IntegrationType I, typename OpBase>
379
380/**
381 * @brief This is specialisation for sources on boundary which depends on normal
382 *
383 * @tparam BASE_DIM
384 * @tparam OpBase
385 */
386template <int FIELD_DIM, IntegrationType I, typename OpBase>
389 : public OpNormalMixVecTimesScalarImpl<FIELD_DIM, I, OpBase> {
391 OpBase>::OpNormalMixVecTimesScalarImpl;
392};
393
394template <typename OpBase>
397 const std::string field_name,
398 ScalarFun source_fun = [](double, double, double) constexpr { return 1; },
399 boost::shared_ptr<Range> ents_ptr = nullptr)
400 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
401 sourceFun(source_fun) {}
402
403protected:
407};
408
409template <typename OpBase>
412 const std::string field_name,
413 ScalarFun source_fun = [](double, double, double) constexpr { return 1; },
414 boost::shared_ptr<Range> ents_ptr = nullptr)
415 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
416 sourceFun(source_fun) {}
417
418protected:
422};
423
424/**
425 * @brief Multiply vector times normal on the face times vector field
426 *
427 * This operator typically will be used to evaluate natural boundary conditions
428 * for mixed formulation.
429 *
430 * @tparam BASE_DIM
431 * @tparam SPACE_DIM
432 * @tparam OpBase
433 */
434template <int SPACE_DIM, IntegrationType I, typename OpBase>
436
437template <int SPACE_DIM, typename OpBase>
439 : public OpBase {
441 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
442 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; },
443 boost::shared_ptr<Range> ents_ptr = nullptr)
444 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr), uPtr(u_ptr),
445 betaCoeff(beta_coeff) {}
446
447protected:
448 boost::shared_ptr<MatrixDouble> uPtr;
451};
452
453template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
454 typename OpBase>
456
457template <int SPACE_DIM, typename OpBase>
460 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
461 boost::shared_ptr<MatrixDouble> y_grad_ptr,
462 ConstantFun source_fun = []() constexpr { return 1; })
463 : OpBase(field_name, field_name, OpBase::OPROW), uPtr(u_ptr),
464 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
465
466protected:
467 boost::shared_ptr<MatrixDouble> uPtr;
468 boost::shared_ptr<MatrixDouble> yGradPtr;
471};
472
473template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
475 : public OpBase {
477 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
478 boost::shared_ptr<MatrixDouble> y_grad_ptr,
479 ConstantFun source_fun = []() constexpr { return 1; })
480 : OpBase(field_name, field_name, OpBase::OPROW), uPtr(u_ptr),
481 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
482
483protected:
484 boost::shared_ptr<MatrixDouble> uPtr;
485 boost::shared_ptr<MatrixDouble> yGradPtr;
488};
489
490template <typename OpBase>
493 EntitiesFieldData::EntData &row_data) {
495
496 // get element volume
497 const double vol = OpBase::getMeasure();
498 // get integration weights
499 auto t_w = OpBase::getFTensor0IntegrationWeight();
500 // get base function gradient on rows
501 auto t_row_base = row_data.getFTensor0N();
502 // get coordinate at integration points
503 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
504 // loop over integration points
505 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
506 // take into account Jacobian
507 const double alpha =
508 t_w * vol * sourceFun(t_coords(0), t_coords(1), t_coords(2));
509 // loop over rows base functions
510 int rr = 0;
511 for (; rr != OpBase::nbRows; ++rr) {
512 OpBase::locF[rr] += alpha * t_row_base;
513 ++t_row_base;
514 }
515 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
516 ++t_row_base;
517 ++t_coords;
518 ++t_w; // move to another integration weight
519 }
521}
522
523template <int FIELD_DIM, typename OpBase>
526 iNtegrate(EntitiesFieldData::EntData &row_data) {
529
530 // get element volume
531 const double vol = OpBase::getMeasure();
532 // get integration weights
533 auto t_w = OpBase::getFTensor0IntegrationWeight();
534 // get base function gradient on rows
535 auto t_row_base = row_data.getFTensor0N();
536 // get coordinate at integration points
537 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
538 // loop over integration points
539 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
540 // source file
541 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
542 // take into account Jacobian
543 const double alpha = t_w * vol;
544 // loop over rows base functions
546 int rr = 0;
547 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
548 t_nf(i) += alpha * t_row_base * t_source(i);
549 ++t_row_base;
550 ++t_nf;
551 }
552 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
553 ++t_row_base;
554 ++t_coords;
555 ++t_w; // move to another integration weight
556 }
558}
559
560template <int FIELD_DIM, typename OpBase>
563 iNtegrate(EntitiesFieldData::EntData &row_data) {
566
567 const size_t nb_base_functions = row_data.getN().size2() / 3;
568 // get element volume
569 const double vol = OpBase::getMeasure();
570 // get integration weights
571 auto t_w = OpBase::getFTensor0IntegrationWeight();
572 // get base function gradient on rows
573 auto t_row_base = row_data.getFTensor1N<3>();
574 // get coordinate at integration points
575 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
576 // loop over integration points
577 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
578 // source file
579 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
580 // take into account Jacobian
581 const double alpha = t_w * vol;
582 // loop over rows base functions
583 int rr = 0;
584 for (; rr != OpBase::nbRows; ++rr) {
585 OpBase::locF[rr] += alpha * t_row_base(i) * t_source(i);
586 ++t_row_base;
587 }
588 for (; rr < nb_base_functions; ++rr)
589 ++t_row_base;
590 ++t_coords;
591 ++t_w; // move to another integration weight
592 }
594}
595
596template <int S, typename OpBase>
598 EntitiesFieldData::EntData &row_data) {
600
601 // get element volume
602 const double vol = OpBase::getMeasure();
603 // get integration weights
604 auto t_w = OpBase::getFTensor0IntegrationWeight();
605 // get base function gradient on rows
606 auto t_row_base = row_data.getFTensor0N();
607 // get vector values
608 auto t_vec = getFTensor0FromVec<S>(*sourceVec);
609
610#ifndef NDEBUG
611 if (sourceVec->size() != OpBase::nbIntegrationPts) {
612 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
613 "Wrong number of integration points %d != %ld",
614 OpBase::nbIntegrationPts, sourceVec->size());
615 }
616#endif
617
618 // get coords
619 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
620 // loop over integration points
621 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
622 // take into account Jacobian
623 const double alpha =
624 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
625 // loop over rows base functions
626 int rr = 0;
627 for (; rr != OpBase::nbRows; ++rr) {
628 OpBase::locF[rr] += alpha * t_row_base * t_vec;
629 ++t_row_base;
630 }
631 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
632 ++t_row_base;
633 ++t_w; // move to another integration weight
634 ++t_vec;
635 ++t_coords;
636 }
638}
639
640template <int FIELD_DIM, int S, typename OpBase>
642 EntitiesFieldData::EntData &row_data) {
644
645 // get element volume
646 const double vol = OpBase::getMeasure();
647 // get integration weights
648 auto t_w = OpBase::getFTensor0IntegrationWeight();
649 // get base function gradient on rows
650 auto t_row_base = row_data.getFTensor0N();
651 // get coords
652 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
653 // get vector values
654 auto t_vec = getFTensor1FromMat<FIELD_DIM, S>(*sourceVec);
655 // loop over integration points
656 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
657 // take into account Jacobian
658 const double alpha =
659 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
660 // get loc vector tensor
661 auto t_nf = OpBase::template getNf<FIELD_DIM>();
662 // loop over rows base functions
663 int rr = 0;
664 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
665 t_nf(i) += alpha * t_row_base * t_vec(i);
666 ++t_row_base;
667 ++t_nf;
668 }
669 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
670 ++t_row_base;
671 ++t_w; // move to another integration weight
672 ++t_vec;
673 ++t_coords;
674 }
676}
677
678template <int FIELD_DIM, int S, typename OpBase>
680 EntitiesFieldData::EntData &row_data) {
682
683 const size_t nb_base_functions = row_data.getN().size2() / 3;
684 // get element volume
685 const double vol = OpBase::getMeasure();
686 // get integration weights
687 auto t_w = OpBase::getFTensor0IntegrationWeight();
688 // get base function gradient on rows
689 auto t_row_base = row_data.getFTensor1N<3>();
690 // get coords
691 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
692 // get vector values
693 auto t_vec = getFTensor1FromMat<FIELD_DIM, S>(*sourceVec);
694 // loop over integration points
695 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
696 // take into account Jacobian
697 const double alpha =
698 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
699 // loop over rows base functions
700 int rr = 0;
701 for (; rr != OpBase::nbRows; ++rr) {
702 OpBase::locF[rr] += alpha * t_row_base(i) * t_vec(i);
703 ++t_row_base;
704 }
705 for (; rr < nb_base_functions; ++rr)
706 ++t_row_base;
707 ++t_w; // move to another integration weight
708 ++t_vec;
709 ++t_coords;
710 }
712}
713
714template <int SPACE_DIM, int S, typename OpBase>
717 EntitiesFieldData::EntData &row_data) {
719
720 // get element volume
721 const double vol = OpBase::getMeasure();
722 // get integration weights
723 auto t_w = OpBase::getFTensor0IntegrationWeight();
724 // get base function gradient on rows
725 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
726 // get field gradient values
727 auto t_val_grad = getFTensor1FromMat<SPACE_DIM, S>(*(matVals));
728 // get coordinate at integration points
729 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
730 // loop over integration points
731 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
732 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
733 // take into account Jacobian
734 const double alpha = t_w * beta;
735 // loop over rows base functions
736 int rr = 0;
737 for (; rr != OpBase::nbRows; rr++) {
738 // calculate element of local matrix
739 OpBase::locF[rr] += alpha * (t_row_grad(i) * t_val_grad(i));
740 ++t_row_grad; // move to another element of gradient of base
741 // function on row
742 }
743 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
744 ++t_row_grad;
745
746 ++t_coords;
747 ++t_val_grad;
748 ++t_w; // move to another integration weight
749 }
751}
752
753template <int SPACE_DIM, int S, typename OpBase>
756 EntitiesFieldData::EntData &row_data) {
758
759 // get element volume
760 const double vol = OpBase::getMeasure();
761 // get integration weights
762 auto t_w = OpBase::getFTensor0IntegrationWeight();
763 // get base function gradient on rows
764 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
765 // get field gradient values
766 auto t_val_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
767 // get coordinate at integration points
768 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
769 // loop over integration points
770 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
771 // take into account Jacobian
772 const double alpha =
773 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
774 // get rhs vector
775 auto t_nf = OpBase::template getNf<SPACE_DIM>();
776 // loop over rows base functions
777 int rr = 0;
778 for (; rr != OpBase::nbRows / SPACE_DIM; rr++) {
779 // calculate element of local matrix
780 t_nf(i) += alpha * (t_row_grad(j) * t_val_grad(i, j));
781 ++t_row_grad; // move to another element of gradient of base
782 // function on row
783 ++t_nf;
784 }
785 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
786 ++t_row_grad;
787
788 ++t_coords;
789 ++t_val_grad;
790 ++t_w; // move to another integration weight
791 }
793}
794
795template <int SPACE_DIM, int S, typename OpBase>
798 EntitiesFieldData::EntData &row_data) {
800 // get element volume
801 const double vol = OpBase::getMeasure();
802 // get coords
803 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
804 // get integration weights
805 auto t_w = OpBase::getFTensor0IntegrationWeight();
806 // get base function gradient on rows
807 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
808 // get field gradient values
809 auto t_val_mat = getFTensor2SymmetricFromMat<SPACE_DIM, S>(*(matVals));
810 // loop over integration points
811 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
812 // take into account Jacobian
813 const double alpha =
814 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
815 // get rhs vector
816 auto t_nf = OpBase::template getNf<SPACE_DIM>();
817 // loop over rows base functions
818 int rr = 0;
819 for (; rr != OpBase::nbRows / SPACE_DIM; rr++) {
820 // calculate element of local matrix
821 t_nf(j) += alpha * (t_row_grad(i) * t_val_mat(i, j));
822 ++t_row_grad; // move to another element of gradient of base
823 // function on row
824 ++t_nf;
825 }
826 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
827 ++t_row_grad;
828 ++t_val_mat;
829 ++t_coords;
830 ++t_w; // move to another integration weight
831 }
833}
834
835template <int FIELD_DIM, int SPACE_DIM, typename OpBase,
836 CoordinateTypes CoordSys>
839 EntitiesFieldData::EntData &row_data) {
841
842 const size_t nb_base_functions = row_data.getN().size2() / 3;
843 auto t_w = this->getFTensor0IntegrationWeight();
844 // get coordinate at integration points
845 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
846 auto t_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
847 auto t_base = row_data.getFTensor1N<3>();
848 auto t_u = getFTensor1FromMat<FIELD_DIM>(*(matVals));
849
850 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
851
852 const double alpha = this->getMeasure() * t_w *
853 betaConst(t_coords(0), t_coords(1), t_coords(2));
854 auto t_nf = OpBase::template getNf<FIELD_DIM>();
855
856 size_t bb = 0;
857 for (; bb != this->nbRows / FIELD_DIM; ++bb) {
858 const double t_div_base = t_diff_base(j, j);
859 t_nf(i) += alpha * t_div_base * t_u(i);
860 if constexpr (CoordSys == CYLINDRICAL) {
861 t_nf(i) += t_base(0) * (alpha / t_coords(0)) * t_u(i);
862 }
863 ++t_nf;
864 ++t_diff_base;
865 ++t_base;
866 }
867 for (; bb < nb_base_functions; ++bb) {
868 ++t_diff_base;
869 ++t_base;
870 }
871
872 ++t_u;
873 ++t_w;
874 ++t_coords;
875 }
876
878}
879
880template <int SPACE_DIM, typename OpBase, CoordinateTypes CoordSys>
883 EntitiesFieldData::EntData &row_data) {
885
886 const size_t nb_base_functions = row_data.getN().size2() / 3;
887 auto t_w = this->getFTensor0IntegrationWeight();
888 // get coordinate at integration points
889 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
890 auto t_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
891 auto t_u = getFTensor0FromVec(*(vecVals));
892
893 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
894
895 const double alpha = this->getMeasure() * t_w *
896 betaConst(t_coords(0), t_coords(1), t_coords(2));
897 ;
898
899 size_t bb = 0;
900 for (; bb != this->nbRows; ++bb) {
901 const double t_div_base = t_diff_base(j, j);
902 OpBase::locF[bb] += alpha * t_div_base * t_u;
903 ++t_diff_base;
904 }
905 for (; bb < nb_base_functions; ++bb)
906 ++t_diff_base;
907
908 ++t_u;
909 ++t_w;
910 ++t_coords;
911 }
912
914}
915
916/**
917 * @brief div U times vector
918 *
919 * \f[
920 * \delta u_j = \phi^m\delta\overline{u}^m_j\\
921 * \delta u_{j,i} = \phi^m_{,i}\delta\overline{u}^m_j\\
922 * \textrm{tr}[\delta u_{j,i}] = \delta u_{j,i}\delta_{ji}\\
923 * (\textrm{tr}[\delta u_{j,i}], v) =\\
924 * (\delta u_{j,i} \delta_{ij}, v) =\\
925 * (\delta u_{j,i}, \delta_{ij} v) =\\
926 * (\phi^m_{,i}\delta\overline{u}^m_j, \delta_{ij} v) \\
927 * f_i^m=(\phi^m_{,i}, v)
928 * \f]
929 *
930 * @tparam FIELD_DIM
931 * @tparam SPACE_DIM
932 * @tparam OpBase
933 * @param row_data
934 * @return MoFEMErrorCode
935 */
936template <int FIELD_DIM, typename OpBase, CoordinateTypes CoordSys>
939 EntitiesFieldData::EntData &row_data) {
942
943 // get element volume
944 const double vol = OpBase::getMeasure();
945 // get integration weights
946 auto t_w = OpBase::getFTensor0IntegrationWeight();
947 // get coordinate at integration points
948 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
949 // get base function gradient on rows
950 auto t_row_grad = row_data.getFTensor1DiffN<FIELD_DIM>();
951 // get vector values
952 auto t_vec = getFTensor0FromVec(*sourceVec);
953 // loop over integration points
954 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
955 // take into account Jacobian
956 const double alpha =
957 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
959 // loop over rows base functions
960 int rr = 0;
961 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
962 t_nf(i) += alpha * t_row_grad(i) * t_vec;
963 ++t_row_grad;
964 ++t_nf;
965 }
966 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
967 ++t_row_grad;
968 ++t_w; // move to another integration weight
969 ++t_vec;
970 ++t_coords;
971 }
973}
974
975template <int SPACE_DIM, typename OpBase>
977 EntitiesFieldData::EntData &row_data) {
979
980 const size_t nb_base_functions = row_data.getN().size2() / 3;
981 auto t_w = this->getFTensor0IntegrationWeight();
982 auto t_coords = this->getFTensor1CoordsAtGaussPts();
983 auto t_base = row_data.getFTensor1N<3>();
984 auto t_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
985
986 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
987
988 const double alpha = this->getMeasure() * t_w;
989 auto t_nf = OpBase::template getNf<SPACE_DIM>();
990
991 size_t bb = 0;
992 for (; bb != this->nbRows / SPACE_DIM; ++bb) {
993 t_nf(i) += alpha * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) *
994 t_base(j) * t_grad(i, j);
995 ++t_nf;
996 ++t_base;
997 }
998 for (; bb < nb_base_functions; ++bb)
999 ++t_base;
1000
1001 ++t_grad;
1002 ++t_coords;
1003 ++t_w;
1004 }
1005
1007}
1008
1009template <int SPACE_DIM, typename OpBase>
1011 EntitiesFieldData::EntData &row_data) {
1013
1014 const size_t nb_base_functions = row_data.getN().size2();
1015 auto t_w = this->getFTensor0IntegrationWeight();
1016 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1017 auto t_base = row_data.getFTensor0N();
1018 auto t_div = getFTensor1FromMat<SPACE_DIM>(*(matVals));
1019
1020 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1021 const double alpha = this->getMeasure() * t_w;
1022 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1023
1024 size_t bb = 0;
1025 for (; bb != this->nbRows / SPACE_DIM; ++bb) {
1026 t_nf(i) += alpha * t_base *
1027 betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * t_div(i);
1028 ++t_nf;
1029 ++t_base;
1030 }
1031 for (; bb < nb_base_functions; ++bb)
1032 ++t_base;
1033
1034 ++t_coords;
1035 ++t_div;
1036 ++t_w;
1037 }
1038
1040}
1041
1042template <typename OpBase>
1044 EntitiesFieldData::EntData &row_data) {
1046
1047 const size_t nb_base_functions = row_data.getN().size2() / 3;
1048 // get element volume
1049 // get integration weights
1050 auto t_w = OpBase::getFTensor0IntegrationWeight();
1051 // get base function gradient on rows
1052 auto t_row_base = row_data.getFTensor1N<3>();
1053 // get coordinate at integration points
1054 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1055 // get normal
1056 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1057 double a = 1;
1058 if (this->getNumeredEntFiniteElementPtr()->getEntType() == MBTRI)
1059 a *= 2;
1060 // loop over integration points
1061 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1062 // take into account Jacobian
1063 const double alpha =
1064 t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2)) / a;
1065 // loop over rows base functions
1066 int rr = 0;
1067 for (; rr != OpBase::nbRows; ++rr) {
1068 OpBase::locF[rr] += alpha * t_row_base(i) * t_normal(i);
1069 ++t_row_base;
1070 }
1071 for (; rr < nb_base_functions; ++rr)
1072 ++t_row_base;
1073 ++t_coords;
1074 ++t_w; // move to another integration weight
1075 ++t_normal;
1076 }
1078}
1079
1080template <typename OpBase>
1082 EntitiesFieldData::EntData &row_data) {
1084
1085 const size_t nb_base_functions = row_data.getN().size2() / 3;
1086 FTensor::Tensor1<double, 3> t_z{0., 0., 1.};
1087 // get element volume
1088 // get integration weights
1089 auto t_w = OpBase::getFTensor0IntegrationWeight();
1090 // get base function gradient on rows
1091 auto t_row_base = row_data.getFTensor1N<3>();
1092 // get coordinate at integration points
1093 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1094 // get normal
1095 auto t_tangent = OpBase::getFTensor1TangentAtGaussPts();
1096 // loop over integration points
1097 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1098 // take into account Jacobian
1099 const double alpha = t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2));
1101 FTensor::Index<'j', 3> j;
1102 FTensor::Index<'k', 3> k;
1103 t_normal(i) = FTensor::levi_civita(i, j, k) * t_tangent(j) * t_z(k);
1104 int rr = 0;
1105 for (; rr != OpBase::nbRows; ++rr) {
1106 OpBase::locF[rr] += alpha * t_row_base(i) * t_normal(i);
1107 ++t_row_base;
1108 }
1109 for (; rr < nb_base_functions; ++rr)
1110 ++t_row_base;
1111 ++t_coords;
1112 ++t_tangent;
1113 ++t_w; // move to another integration weight
1114 }
1116}
1117
1118template <int SPACE_DIM, typename OpBase>
1121 EntitiesFieldData::EntData &row_data) {
1123 FTensor::Index<'i', 3> i;
1125
1126 const size_t nb_base_functions = row_data.getN().size2() / 3;
1127 // get element volume
1128 // get integration weights
1129 auto t_w = OpBase::getFTensor0IntegrationWeight();
1130 // get base function gradient on rows
1131 auto t_row_base = row_data.getFTensor1N<3>();
1132 // get coordinate at integration points
1133 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1134 // get normal
1135 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1136 // get field
1137 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1138 // loop over integration points
1139 auto a = OpBase::getMeasure();
1140 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1141 // take into account Jacobian
1142 auto l2 = std::sqrt(t_normal(i) * t_normal(i));
1143 const double alpha =
1144 t_w * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * (a / l2);
1145 // get rhs vector
1146 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1147 // loop over rows base functions
1148 int rr = 0;
1149 for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1150 t_nf(J) += alpha * (t_row_base(i) * t_normal(i)) * t_u(J);
1151 ++t_row_base;
1152 ++t_nf;
1153 }
1154 for (; rr < nb_base_functions; ++rr)
1155 ++t_row_base;
1156 ++t_coords;
1157 ++t_w; // move to another integration weight
1158 ++t_normal;
1159 ++t_u;
1160 }
1162}
1163
1164template <int SPACE_DIM, typename OpBase>
1167 EntitiesFieldData::EntData &row_data) {
1169
1170 auto t_w = this->getFTensor0IntegrationWeight();
1171 auto t_base = row_data.getFTensor0N();
1172
1173 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1174 auto t_grad_y = getFTensor1FromMat<SPACE_DIM>(*yGradPtr);
1175
1177 const double alpha_constant = alphaConstant();
1178 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1179
1180 // get element volume
1181 const double vol = OpBase::getMeasure();
1182 const double c = (t_grad_y(i) * t_u(i)) * (t_w * vol * alpha_constant);
1183
1184 // get element volume
1185 int rr = 0;
1186 for (; rr != OpBase::nbRows; ++rr) {
1187 OpBase::locF[rr] += c * t_base;
1188 ++t_base;
1189 }
1190 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1191 ++t_base;
1192
1193 ++t_w; // move to another integration weight
1194 ++t_u;
1195 ++t_grad_y;
1196 }
1197
1199}
1200
1201template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
1204 EntitiesFieldData::EntData &row_data) {
1206
1207 auto t_w = this->getFTensor0IntegrationWeight();
1208 auto t_base = row_data.getFTensor0N();
1209
1210 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1211 auto t_grad_y = getFTensor2FromMat<FIELD_DIM, SPACE_DIM>(*yGradPtr);
1212
1215 const double alpha_constant = alphaConstant();
1216 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1217
1218 // get element volume
1219 const double vol = OpBase::getMeasure();
1220
1222 t_c(J) = (t_grad_y(J, i) * t_u(i)) * (t_w * vol * alpha_constant);
1223
1225 int rr = 0;
1226 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
1227 t_nf(J) += t_base * t_c(J);
1228 ++t_base;
1229 ++t_nf;
1230 }
1231 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1232 ++t_base;
1233 ++t_w; // move to another integration weight
1234
1235 ++t_u;
1236 ++t_grad_y;
1237 }
1238
1240}
1241
1242} // namespace MoFEM
1243
1244#endif // __LINEAR_FORMS_INTEGRATORS_HPP__
constexpr double a
constexpr int SPACE_DIM
constexpr int FIELD_DIM
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
CoordinateTypes
Coodinate system.
@ CYLINDRICAL
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
constexpr int BASE_DIM
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.
FTensor::Index< 'i', SPACE_DIM > i
const double c
speed of light (cm/ns)
FTensor::Index< 'J', DIM1 > J
Definition level_set.cpp:30
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.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
boost::function< double(double)> TimeFun
Lambda function used to scale with time.
FTensor::Tensor1< FTensor::PackPtr< T *, S >, Tensor_Dim > getFTensor1FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 1 (vector) form data matrix.
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, Tensor_Dim1, Tensor_Dim2 > getFTensor2FromMat(MatrixDouble &data)
Get tensor rank 2 (matrix) form data matrix.
static auto getFTensor2SymmetricFromMat(ublas::matrix< T, L, A > &data)
Get symmetric tensor rank 2 (matrix) form data matrix.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
boost::function< double()> ConstantFun
Constant function type.
auto getFTensor1FromArray(VectorDouble &data)
Get FTensor1 from array.
constexpr IntegrationType I
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
OpBaseTimesScalarImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpBaseTimesVectorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > vec, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, 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=[]() constexpr { return 1;})
OpConvectiveTermRhsImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, boost::shared_ptr< MatrixDouble > y_grad_ptr, ConstantFun source_fun=[]() constexpr { return 1;})
OpGradTimesSymTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;})
OpGradTimesTensorImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, 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) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
MoFEMErrorCode iNtegrate(EntitiesFieldData::EntData &data)
Class dedicated to integrate operator.
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec, ScalarFun beta=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< VectorDouble > vec_vals, ScalarFun beta=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixDivTimesUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta=[](double, double, double) { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpMixTensorTimesGradUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals)
OpMixTensorTimesGradUImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_fun)
Tensor field time gradient of vector field.
OpMixVecTimesDivLambdaImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals, ScalarFun beta_fun)
OpMixVecTimesDivLambdaImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > mat_vals)
Vector field time divergence of tensor.
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
OpNormalMixVecTimesScalarImpl(const std::string field_name, ScalarFun source_fun=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
Multiply vector times normal on the face times scalar function.
OpNormalMixVecTimesVectorFieldImpl(const std::string field_name, boost::shared_ptr< MatrixDouble > u_ptr, ScalarFun beta_coeff=[](double, double, double) constexpr { return 1;}, boost::shared_ptr< Range > ents_ptr=nullptr)
Multiply vector times normal on the face times vector field.
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, TimeFun time_fun, 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)
Construct a new Op Source Impl object.
OpSourceImpl(const std::string field_name, TimeFun time_fun, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
Construct a new Op Source Impl object.
OpSourceImpl(const std::string field_name, TimeFun time_fun, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)
OpSourceImpl(const std::string field_name, VectorFun< FIELD_DIM > source_fun, boost::shared_ptr< Range > ents_ptr=nullptr)