v0.14.0
Loading...
Searching...
No Matches
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#ifndef __LINEAR_FORMS_INTEGRATORS_HPP__
8#define __LINEAR_FORMS_INTEGRATORS_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
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
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
298template <int SPACE_DIM, IntegrationType I, typename OpBase>
300
301template <int SPACE_DIM, typename OpBase>
304 boost::shared_ptr<MatrixDouble> mat_vals)
305 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals) {}
306
308 boost::shared_ptr<MatrixDouble> mat_vals,
309 ScalarFun beta_fun)
310 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
311 betaCoeff(beta_fun) {}
312
313protected:
314 boost::shared_ptr<MatrixDouble> matVals;
315 ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
316
319};
320
321template <int SPACE_DIM, IntegrationType I, typename OpBase>
323
324template <int SPACE_DIM, typename OpBase>
327 boost::shared_ptr<MatrixDouble> mat_vals)
328 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals) {}
329
331 boost::shared_ptr<MatrixDouble> mat_vals,
332 ScalarFun beta_fun)
333 : OpBase(field_name, field_name, OpBase::OPROW), matVals(mat_vals),
334 betaCoeff(beta_fun) {}
335
336protected:
337 boost::shared_ptr<MatrixDouble> matVals;
338 ScalarFun betaCoeff = [](double, double, double) constexpr { return 1; };
339
343};
344
345/**
346 * @brief Multiply vector times normal on the face times scalar function
347 *
348 * This operator typically will be used to evaluate natural boundary conditions
349 * for mixed formulation.
350 *
351 * @tparam BASE_DIM
352 * @tparam SPACE_DIM
353 * @tparam OpBase
354 */
355template <int SPACE_DIM, IntegrationType I, typename OpBase>
357
358/**
359 * @brief This is specialisation for sources on boundary which depends on normal
360 *
361 * @tparam BASE_DIM
362 * @tparam OpBase
363 */
364template <int FIELD_DIM, IntegrationType I, typename OpBase>
367 : public OpNormalMixVecTimesScalarImpl<FIELD_DIM, I, OpBase> {
369 OpBase>::OpNormalMixVecTimesScalarImpl;
370};
371
372template <typename OpBase>
375 const std::string field_name,
376 ScalarFun source_fun = [](double, double, double) constexpr { return 1; },
377 boost::shared_ptr<Range> ents_ptr = nullptr)
378 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
379 sourceFun(source_fun) {}
380
381protected:
385};
386
387template <typename OpBase>
390 const std::string field_name,
391 ScalarFun source_fun = [](double, double, double) constexpr { return 1; },
392 boost::shared_ptr<Range> ents_ptr = nullptr)
393 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr),
394 sourceFun(source_fun) {}
395
396protected:
400};
401
402/**
403 * @brief Multiply vector times normal on the face times vector field
404 *
405 * This operator typically will be used to evaluate natural boundary conditions
406 * for mixed formulation.
407 *
408 * @tparam BASE_DIM
409 * @tparam SPACE_DIM
410 * @tparam OpBase
411 */
412template <int SPACE_DIM, IntegrationType I, typename OpBase>
414
415template <int SPACE_DIM, typename OpBase>
417 : public OpBase {
419 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
420 ScalarFun beta_coeff = [](double, double, double) constexpr { return 1; },
421 boost::shared_ptr<Range> ents_ptr = nullptr)
422 : OpBase(field_name, field_name, OpBase::OPROW, ents_ptr), uPtr(u_ptr),
423 betaCoeff(beta_coeff) {}
424
425protected:
426 boost::shared_ptr<MatrixDouble> uPtr;
429};
430
431template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, IntegrationType I,
432 typename OpBase>
434
435template <int SPACE_DIM, typename OpBase>
438 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
439 boost::shared_ptr<MatrixDouble> y_grad_ptr,
440 ConstantFun source_fun = []() constexpr { return 1; })
441 : OpBase(field_name, field_name, OpBase::OPROW), uPtr(u_ptr),
442 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
443
444protected:
445 boost::shared_ptr<MatrixDouble> uPtr;
446 boost::shared_ptr<MatrixDouble> yGradPtr;
449};
450
451template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
453 : public OpBase {
455 const std::string field_name, boost::shared_ptr<MatrixDouble> u_ptr,
456 boost::shared_ptr<MatrixDouble> y_grad_ptr,
457 ConstantFun source_fun = []() constexpr { return 1; })
458 : OpBase(field_name, field_name, OpBase::OPROW), uPtr(u_ptr),
459 yGradPtr(y_grad_ptr), alphaConstant(source_fun) {}
460
461protected:
462 boost::shared_ptr<MatrixDouble> uPtr;
463 boost::shared_ptr<MatrixDouble> yGradPtr;
466};
467
468/**
469 * @brief Linear integrator form
470 * @ingroup mofem_forms
471 *
472 * @tparam EleOp
473 * @tparam A
474 * @tparam I
475 */
476template <typename EleOp>
477template <AssemblyType A>
478template <IntegrationType I>
479struct FormsIntegrators<EleOp>::Assembly<A>::LinearForm {
480
481 /**
482 * @brief Integrate \f$(v,f(\mathbf{x}))_\Omega\f$, f is a scalar
483 * @ingroup mofem_forms
484 *
485 * @note \f$f(\mathbf{x})\$ is scalar function of coordinates
486 *
487 * @tparam BASE_DIM
488 * @tparam FIELD_DIM
489 */
490 template <int BASE_DIM, int FIELD_DIM,
491 typename S = SourceFunctionSpecialization>
492 using OpSource =
494
495 /**
496 * @brief Vector field integrator \f$(v_i,f)_\Omega\f$, f is a vector
497 * @ingroup mofem_forms
498 *
499 * @tparam BASE_DIM
500 */
501 template <int BASE_DIM, int S = 1>
502 using OpBaseTimesScalar = OpBaseTimesScalarImpl<BASE_DIM, S, I, OpBase>;
503
504 /** @deprecated use instead OpBaseTimesScalar
505 */
506 template <int BASE_DIM, int S = 1>
507 using OpBaseTimesScalarField = OpBaseTimesScalar<BASE_DIM, S>;
508
509 /**
510 * @brief Vector field integrator \f$(v,f_i)_\Omega\f$, f is a vector
511 * @ingroup mofem_forms
512 *
513 * @tparam BASE_DIM
514 * @tparam FIELD_DIM
515 * @tparam 0
516 */
517 template <int BASE_DIM, int FIELD_DIM, int S>
518 using OpBaseTimesVector =
520 //! [Source operator]
521
522 //! [Grad times tensor]
523
524 /**
525 * @brief Integrate \f$(v_{,i},f_{ij})_\Omega\f$, f is a vector
526 * @ingroup mofem_forms
527 *
528 * @note \f$f_{ij}\$ is tensor at integration points
529 *
530 * @tparam BASE_DIM
531 * @tparam FIELD_DIM
532 * @tparam SPACE_DIM
533 */
534 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
537
538 /**
539 * @brief Integrate \f$(v_{,i},f_{ij})_\Omega\f$, f is symmetric tensor
540 * @ingroup mofem_forms
541 *
542 * @note \f$f_{ij}\$ is tensor at integration points
543 *
544 * @tparam BASE_DIM
545 * @tparam FIELD_DIM
546 * @tparam SPACE_DIM
547 */
548 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM, int S = 1>
549 using OpGradTimesSymTensor =
551
552 /**
553 * @brief Integrate \f$(\lambda_{ij,j},u_{i})_\Omega\f$
554 *
555 * @tparam BASE_DIM
556 * @tparam FIELD_DIM
557 * @tparam SPACE_DIM
558 */
559 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM,
560 CoordinateTypes CoordSys = CARTESIAN>
561 using OpMixDivTimesU =
563
564 /**
565 * @brief Integrate \f$(\lambda_{ij},u_{i,j})_\Omega\f$
566 *
567 * @tparam SPACE_DIM
568 */
569 template <int SPACE_DIM>
570 using OpMixTensorTimesGradU = OpMixTensorTimesGradUImpl<SPACE_DIM, I, OpBase>;
571
572 /**
573 * @brief Integrate \f$(u_{i},\lambda_{ij,j})_\Omega\f$
574 *
575 * @tparam SPACE_DIM
576 */
577 template <int SPACE_DIM>
578 using OpMixVecTimesDivLambda =
580
581 /**
582 * @brief Multiply vector times normal on the face times scalar function
583 *
584 * This operator typically will be used to evaluate natural boundary
585 * conditions for mixed formulation.
586 *
587 * @tparam BASE_DIM
588 * @tparam SPACE_DIM
589 * @tparam OpBase
590 */
591 template <int SPACE_DIM>
592 using OpNormalMixVecTimesScalar =
594
595 /**
596 * @brief Multiply vector times normal on the face times vector field
597 *
598 * This operator typically will be used to evaluate natural boundary
599 * conditions for mixed formulation.
600 *
601 * @tparam BASE_DIM
602 * @tparam SPACE_DIM
603 * @tparam OpBase
604 */
605 template <int SPACE_DIM>
606 using OpNormalMixVecTimesVectorField =
608
609 /**
610 * @brief Convective term
611 *
612 * \f[
613 * (v, u_i \mathbf{y}_{,i})
614 * \f]
615 * where \f$\mathbf{y}\$ can be scalar or vector.
616 *
617 * @tparam BASE_DIM
618 * @tparam FIELD_DIM
619 * @tparam SPACE_DIM
620 */
621 template <int BASE_DIM, int FIELD_DIM, int SPACE_DIM>
624};
625
626template <typename OpBase>
629 EntitiesFieldData::EntData &row_data) {
631
632 // get element volume
633 const double vol = OpBase::getMeasure();
634 // get integration weights
635 auto t_w = OpBase::getFTensor0IntegrationWeight();
636 // get base function gradient on rows
637 auto t_row_base = row_data.getFTensor0N();
638 // get coordinate at integration points
639 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
640 // loop over integration points
641 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
642 // take into account Jacobian
643 const double alpha =
644 t_w * vol * sourceFun(t_coords(0), t_coords(1), t_coords(2));
645 // loop over rows base functions
646 int rr = 0;
647 for (; rr != OpBase::nbRows; ++rr) {
648 OpBase::locF[rr] += alpha * t_row_base;
649 ++t_row_base;
650 }
651 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
652 ++t_row_base;
653 ++t_coords;
654 ++t_w; // move to another integration weight
655 }
657}
658
659template <int FIELD_DIM, typename OpBase>
662 iNtegrate(EntitiesFieldData::EntData &row_data) {
665
666 // get element volume
667 const double vol = OpBase::getMeasure();
668 // get integration weights
669 auto t_w = OpBase::getFTensor0IntegrationWeight();
670 // get base function gradient on rows
671 auto t_row_base = row_data.getFTensor0N();
672 // get coordinate at integration points
673 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
674 // loop over integration points
675 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
676 // source file
677 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
678 // take into account Jacobian
679 const double alpha = t_w * vol;
680 // loop over rows base functions
681 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(OpBase::locF);
682 int rr = 0;
683 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
684 t_nf(i) += alpha * t_row_base * t_source(i);
685 ++t_row_base;
686 ++t_nf;
687 }
688 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
689 ++t_row_base;
690 ++t_coords;
691 ++t_w; // move to another integration weight
692 }
694}
695
696template <int FIELD_DIM, typename OpBase>
699 iNtegrate(EntitiesFieldData::EntData &row_data) {
702
703 const size_t nb_base_functions = row_data.getN().size2() / 3;
704 // get element volume
705 const double vol = OpBase::getMeasure();
706 // get integration weights
707 auto t_w = OpBase::getFTensor0IntegrationWeight();
708 // get base function gradient on rows
709 auto t_row_base = row_data.getFTensor1N<3>();
710 // get coordinate at integration points
711 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
712 // loop over integration points
713 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
714 // source file
715 auto t_source = sourceFun(t_coords(0), t_coords(1), t_coords(2));
716 // take into account Jacobian
717 const double alpha = t_w * vol;
718 // loop over rows base functions
719 int rr = 0;
720 for (; rr != OpBase::nbRows; ++rr) {
721 OpBase::locF[rr] += alpha * t_row_base(i) * t_source(i);
722 ++t_row_base;
723 }
724 for (; rr < nb_base_functions; ++rr)
725 ++t_row_base;
726 ++t_coords;
727 ++t_w; // move to another integration weight
728 }
730}
731
732template <int S, typename OpBase>
734 EntitiesFieldData::EntData &row_data) {
736
737 // get element volume
738 const double vol = OpBase::getMeasure();
739 // get integration weights
740 auto t_w = OpBase::getFTensor0IntegrationWeight();
741 // get base function gradient on rows
742 auto t_row_base = row_data.getFTensor0N();
743 // get vector values
744 auto t_vec = getFTensor0FromVec<S>(*sourceVec);
745
746#ifndef NDEBUG
747 if (sourceVec->size() != OpBase::nbIntegrationPts) {
748 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
749 "Wrong number of integration points %d != %d",
750 OpBase::nbIntegrationPts, sourceVec->size());
751 }
752#endif
753
754 // get coords
755 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
756 // loop over integration points
757 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
758 // take into account Jacobian
759 const double alpha =
760 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
761 // loop over rows base functions
762 int rr = 0;
763 for (; rr != OpBase::nbRows; ++rr) {
764 OpBase::locF[rr] += alpha * t_row_base * t_vec;
765 ++t_row_base;
766 }
767 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
768 ++t_row_base;
769 ++t_w; // move to another integration weight
770 ++t_vec;
771 ++t_coords;
772 }
774}
775
776template <int FIELD_DIM, int S, typename OpBase>
778 EntitiesFieldData::EntData &row_data) {
780
781 // get element volume
782 const double vol = OpBase::getMeasure();
783 // get integration weights
784 auto t_w = OpBase::getFTensor0IntegrationWeight();
785 // get base function gradient on rows
786 auto t_row_base = row_data.getFTensor0N();
787 // get coords
788 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
789 // get vector values
790 auto t_vec = getFTensor1FromMat<FIELD_DIM, S>(*sourceVec);
791 // loop over integration points
792 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
793 // take into account Jacobian
794 const double alpha =
795 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
796 // get loc vector tensor
797 auto t_nf = OpBase::template getNf<FIELD_DIM>();
798 // loop over rows base functions
799 int rr = 0;
800 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
801 t_nf(i) += alpha * t_row_base * t_vec(i);
802 ++t_row_base;
803 ++t_nf;
804 }
805 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
806 ++t_row_base;
807 ++t_w; // move to another integration weight
808 ++t_vec;
809 ++t_coords;
810 }
812}
813
814template <int FIELD_DIM, int S, typename OpBase>
816 EntitiesFieldData::EntData &row_data) {
818
819 const size_t nb_base_functions = row_data.getN().size2() / 3;
820 // get element volume
821 const double vol = OpBase::getMeasure();
822 // get integration weights
823 auto t_w = OpBase::getFTensor0IntegrationWeight();
824 // get base function gradient on rows
825 auto t_row_base = row_data.getFTensor1N<3>();
826 // get coords
827 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
828 // get vector values
829 auto t_vec = getFTensor1FromMat<FIELD_DIM, S>(*sourceVec);
830 // loop over integration points
831 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
832 // take into account Jacobian
833 const double alpha =
834 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
835 // loop over rows base functions
836 int rr = 0;
837 for (; rr != OpBase::nbRows; ++rr) {
838 OpBase::locF[rr] += alpha * t_row_base(i) * t_vec(i);
839 ++t_row_base;
840 }
841 for (; rr < nb_base_functions; ++rr)
842 ++t_row_base;
843 ++t_w; // move to another integration weight
844 ++t_vec;
845 ++t_coords;
846 }
848}
849
850template <int SPACE_DIM, int S, typename OpBase>
853 EntitiesFieldData::EntData &row_data) {
855
856 // get element volume
857 const double vol = OpBase::getMeasure();
858 // get integration weights
859 auto t_w = OpBase::getFTensor0IntegrationWeight();
860 // get base function gradient on rows
861 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
862 // get filed gradient values
863 auto t_val_grad = getFTensor1FromMat<SPACE_DIM, S>(*(matVals));
864 // get coordinate at integration points
865 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
866 // loop over integration points
867 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
868 const double beta = vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
869 // take into account Jacobian
870 const double alpha = t_w * beta;
871 // loop over rows base functions
872 int rr = 0;
873 for (; rr != OpBase::nbRows; rr++) {
874 // calculate element of local matrix
875 OpBase::locF[rr] += alpha * (t_row_grad(i) * t_val_grad(i));
876 ++t_row_grad; // move to another element of gradient of base
877 // function on row
878 }
879 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
880 ++t_row_grad;
881
882 ++t_coords;
883 ++t_val_grad;
884 ++t_w; // move to another integration weight
885 }
887}
888
889template <int SPACE_DIM, int S, typename OpBase>
892 EntitiesFieldData::EntData &row_data) {
894
895 // get element volume
896 const double vol = OpBase::getMeasure();
897 // get integration weights
898 auto t_w = OpBase::getFTensor0IntegrationWeight();
899 // get base function gradient on rows
900 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
901 // get filed gradient values
902 auto t_val_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
903 // get coordinate at integration points
904 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
905 // loop over integration points
906 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
907 // take into account Jacobian
908 const double alpha =
909 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
910 // get rhs vector
911 auto t_nf = OpBase::template getNf<SPACE_DIM>();
912 // loop over rows base functions
913 int rr = 0;
914 for (; rr != OpBase::nbRows / SPACE_DIM; rr++) {
915 // calculate element of local matrix
916 t_nf(i) += alpha * (t_row_grad(j) * t_val_grad(i, j));
917 ++t_row_grad; // move to another element of gradient of base
918 // function on row
919 ++t_nf;
920 }
921 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
922 ++t_row_grad;
923
924 ++t_coords;
925 ++t_val_grad;
926 ++t_w; // move to another integration weight
927 }
929}
930
931template <int SPACE_DIM, int S, typename OpBase>
934 EntitiesFieldData::EntData &row_data) {
936 // get element volume
937 const double vol = OpBase::getMeasure();
938 // get coords
939 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
940 // get integration weights
941 auto t_w = OpBase::getFTensor0IntegrationWeight();
942 // get base function gradient on rows
943 auto t_row_grad = row_data.getFTensor1DiffN<SPACE_DIM>();
944 // get filed gradient values
945 auto t_val_mat = getFTensor2SymmetricFromMat<SPACE_DIM, S>(*(matVals));
946 // loop over integration points
947 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
948 // take into account Jacobian
949 const double alpha =
950 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
951 // get rhs vector
952 auto t_nf = OpBase::template getNf<SPACE_DIM>();
953 // loop over rows base functions
954 int rr = 0;
955 for (; rr != OpBase::nbRows / SPACE_DIM; rr++) {
956 // calculate element of local matrix
957 t_nf(j) += alpha * (t_row_grad(i) * t_val_mat(i, j));
958 ++t_row_grad; // move to another element of gradient of base
959 // function on row
960 ++t_nf;
961 }
962 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
963 ++t_row_grad;
964 ++t_val_mat;
965 ++t_coords;
966 ++t_w; // move to another integration weight
967 }
969}
970
971template <int FIELD_DIM, int SPACE_DIM, typename OpBase,
972 CoordinateTypes CoordSys>
975 EntitiesFieldData::EntData &row_data) {
977
978 const size_t nb_base_functions = row_data.getN().size2() / 3;
979 auto t_w = this->getFTensor0IntegrationWeight();
980 // get coordinate at integration points
981 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
982 auto t_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
983 auto t_base = row_data.getFTensor1N<3>();
984 auto t_u = getFTensor1FromMat<FIELD_DIM>(*(matVals));
985
986 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
987
988 const double alpha = this->getMeasure() * t_w *
989 betaConst(t_coords(0), t_coords(1), t_coords(2));
990 auto t_nf = OpBase::template getNf<FIELD_DIM>();
991
992 size_t bb = 0;
993 for (; bb != this->nbRows / FIELD_DIM; ++bb) {
994 const double t_div_base = t_diff_base(j, j);
995 t_nf(i) += alpha * t_div_base * t_u(i);
996 if constexpr (CoordSys == CYLINDRICAL) {
997 t_nf(i) += alpha * (t_base(0) / t_coords(0)) * t_u(i);
998 }
999 ++t_nf;
1000 ++t_diff_base;
1001 ++t_base;
1002 }
1003 for (; bb < nb_base_functions; ++bb) {
1004 ++t_diff_base;
1005 ++t_base;
1006 }
1007
1008 ++t_u;
1009 ++t_w;
1010 ++t_coords;
1011 }
1012
1014}
1015
1016template <int SPACE_DIM, typename OpBase, CoordinateTypes CoordSys>
1019 EntitiesFieldData::EntData &row_data) {
1021
1022 const size_t nb_base_functions = row_data.getN().size2() / 3;
1023 auto t_w = this->getFTensor0IntegrationWeight();
1024 // get coordinate at integration points
1025 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1026 auto t_diff_base = row_data.getFTensor2DiffN<3, SPACE_DIM>();
1027 auto t_u = getFTensor0FromVec(*(vecVals));
1028
1029 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1030
1031 const double alpha = this->getMeasure() * t_w *
1032 betaConst(t_coords(0), t_coords(1), t_coords(2));
1033 ;
1034
1035 size_t bb = 0;
1036 for (; bb != this->nbRows; ++bb) {
1037 const double t_div_base = t_diff_base(j, j);
1038 OpBase::locF[bb] += alpha * t_div_base * t_u;
1039 ++t_diff_base;
1040 }
1041 for (; bb < nb_base_functions; ++bb)
1042 ++t_diff_base;
1043
1044 ++t_u;
1045 ++t_w;
1046 ++t_coords;
1047 }
1048
1050}
1051
1052/**
1053 * @brief div U times vector
1054 *
1055 * \f[
1056 * \delta u_j = \phi^m\delta\overline{u}^m_j\\
1057 * \delta u_{j,i} = \phi^m_{,i}\delta\overline{u}^m_j\\
1058 * \textrm{tr}[\delta u_{j,i}] = \delta u_{j,i}\delta_{ji}\\
1059 * (\textrm{tr}[\delta u_{j,i}], v) =\\
1060 * (\delta u_{j,i} \delta_{ij}, v) =\\
1061 * (\delta u_{j,i}, \delta_{ij} v) =\\
1062 * (\phi^m_{,i}\delta\overline{u}^m_j, \delta_{ij} v) \\
1063 * f_i^m=(\phi^m_{,i}, v)
1064 * \f]
1065 *
1066 * @tparam FIELD_DIM
1067 * @tparam SPACE_DIM
1068 * @tparam OpBase
1069 * @param row_data
1070 * @return MoFEMErrorCode
1071 */
1072template <int FIELD_DIM, typename OpBase, CoordinateTypes CoordSys>
1075 EntitiesFieldData::EntData &row_data) {
1078
1079 // get element volume
1080 const double vol = OpBase::getMeasure();
1081 // get integration weights
1082 auto t_w = OpBase::getFTensor0IntegrationWeight();
1083 // get coordinate at integration points
1084 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1085 // get base function gradient on rows
1086 auto t_row_grad = row_data.getFTensor1DiffN<FIELD_DIM>();
1087 // get vector values
1088 auto t_vec = getFTensor0FromVec(*sourceVec);
1089 // loop over integration points
1090 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1091 // take into account Jacobian
1092 const double alpha =
1093 t_w * vol * betaCoeff(t_coords(0), t_coords(1), t_coords(2));
1094 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(OpBase::locF);
1095 // loop over rows base functions
1096 int rr = 0;
1097 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
1098 t_nf(i) += alpha * t_row_grad(i) * t_vec;
1099 ++t_row_grad;
1100 ++t_nf;
1101 }
1102 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1103 ++t_row_grad;
1104 ++t_w; // move to another integration weight
1105 ++t_vec;
1106 ++t_coords;
1107 }
1109}
1110
1111template <int SPACE_DIM, typename OpBase>
1113 EntitiesFieldData::EntData &row_data) {
1115
1116 const size_t nb_base_functions = row_data.getN().size2() / 3;
1117 auto t_w = this->getFTensor0IntegrationWeight();
1118 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1119 auto t_base = row_data.getFTensor1N<3>();
1120 auto t_grad = getFTensor2FromMat<SPACE_DIM, SPACE_DIM>(*(matVals));
1121
1122 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1123
1124 const double alpha = this->getMeasure() * t_w;
1125 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1126
1127 size_t bb = 0;
1128 for (; bb != this->nbRows / SPACE_DIM; ++bb) {
1129 t_nf(i) += alpha * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) *
1130 t_base(j) * t_grad(i, j);
1131 ++t_nf;
1132 ++t_base;
1133 }
1134 for (; bb < nb_base_functions; ++bb)
1135 ++t_base;
1136
1137 ++t_grad;
1138 ++t_coords;
1139 ++t_w;
1140 }
1141
1143}
1144
1145template <int SPACE_DIM, typename OpBase>
1147 EntitiesFieldData::EntData &row_data) {
1149
1150 const size_t nb_base_functions = row_data.getN().size2();
1151 auto t_w = this->getFTensor0IntegrationWeight();
1152 auto t_coords = this->getFTensor1CoordsAtGaussPts();
1153 auto t_base = row_data.getFTensor0N();
1154 auto t_div = getFTensor1FromMat<SPACE_DIM>(*(matVals));
1155
1156 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1157 const double alpha = this->getMeasure() * t_w;
1158 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1159
1160 size_t bb = 0;
1161 for (; bb != this->nbRows / SPACE_DIM; ++bb) {
1162 t_nf(i) += alpha * t_base *
1163 betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * t_div(i);
1164 ++t_nf;
1165 ++t_base;
1166 }
1167 for (; bb < nb_base_functions; ++bb)
1168 ++t_base;
1169
1170 ++t_coords;
1171 ++t_div;
1172 ++t_w;
1173 }
1174
1176}
1177
1178template <typename OpBase>
1180 EntitiesFieldData::EntData &row_data) {
1182
1183 const size_t nb_base_functions = row_data.getN().size2() / 3;
1184 // get element volume
1185 // get integration weights
1186 auto t_w = OpBase::getFTensor0IntegrationWeight();
1187 // get base function gradient on rows
1188 auto t_row_base = row_data.getFTensor1N<3>();
1189 // get coordinate at integration points
1190 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1191 // get normal
1192 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1193 double a = 1;
1194 if (this->getNumeredEntFiniteElementPtr()->getEntType() == MBTRI)
1195 a *= 2;
1196 // loop over integration points
1197 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1198 // take into account Jacobian
1199 const double alpha =
1200 t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2)) / a;
1201 // loop over rows base functions
1202 int rr = 0;
1203 for (; rr != OpBase::nbRows; ++rr) {
1204 OpBase::locF[rr] += alpha * t_row_base(i) * t_normal(i);
1205 ++t_row_base;
1206 }
1207 for (; rr < nb_base_functions; ++rr)
1208 ++t_row_base;
1209 ++t_coords;
1210 ++t_w; // move to another integration weight
1211 ++t_normal;
1212 }
1214}
1215
1216template <typename OpBase>
1218 EntitiesFieldData::EntData &row_data) {
1220
1221 const size_t nb_base_functions = row_data.getN().size2() / 3;
1222 FTensor::Tensor1<double, 3> t_z{0., 0., 1.};
1223 // get element volume
1224 // get integration weights
1225 auto t_w = OpBase::getFTensor0IntegrationWeight();
1226 // get base function gradient on rows
1227 auto t_row_base = row_data.getFTensor1N<3>();
1228 // get coordinate at integration points
1229 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1230 // get normal
1231 auto t_tangent = OpBase::getFTensor1TangentAtGaussPts();
1232 // loop over integration points
1233 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1234 // take into account Jacobian
1235 const double alpha = t_w * sourceFun(t_coords(0), t_coords(1), t_coords(2));
1239 t_normal(i) = FTensor::levi_civita(i, j, k) * t_tangent(j) * t_z(k);
1240 int rr = 0;
1241 for (; rr != OpBase::nbRows; ++rr) {
1242 OpBase::locF[rr] += alpha * t_row_base(i) * t_normal(i);
1243 ++t_row_base;
1244 }
1245 for (; rr < nb_base_functions; ++rr)
1246 ++t_row_base;
1247 ++t_coords;
1248 ++t_tangent;
1249 ++t_w; // move to another integration weight
1250 }
1252}
1253
1254template <int SPACE_DIM, typename OpBase>
1257 EntitiesFieldData::EntData &row_data) {
1261
1262 const size_t nb_base_functions = row_data.getN().size2() / 3;
1263 // get element volume
1264 // get integration weights
1265 auto t_w = OpBase::getFTensor0IntegrationWeight();
1266 // get base function gradient on rows
1267 auto t_row_base = row_data.getFTensor1N<3>();
1268 // get coordinate at integration points
1269 auto t_coords = OpBase::getFTensor1CoordsAtGaussPts();
1270 // get normal
1271 auto t_normal = OpBase::getFTensor1NormalsAtGaussPts();
1272 // get field
1273 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1274 // loop over integration points
1275 auto a = OpBase::getMeasure();
1276 for (int gg = 0; gg != OpBase::nbIntegrationPts; gg++) {
1277 // take into account Jacobian
1278 auto l2 = std::sqrt(t_normal(i) * t_normal(i));
1279 const double alpha =
1280 t_w * betaCoeff(t_coords(0), t_coords(1), t_coords(2)) * (a / l2);
1281 // get rhs vector
1282 auto t_nf = OpBase::template getNf<SPACE_DIM>();
1283 // loop over rows base functions
1284 int rr = 0;
1285 for (; rr != OpBase::nbRows / SPACE_DIM; ++rr) {
1286 t_nf(J) += alpha * (t_row_base(i) * t_normal(i)) * t_u(J);
1287 ++t_row_base;
1288 ++t_nf;
1289 }
1290 for (; rr < nb_base_functions; ++rr)
1291 ++t_row_base;
1292 ++t_coords;
1293 ++t_w; // move to another integration weight
1294 ++t_normal;
1295 ++t_u;
1296 }
1298}
1299
1300template <int SPACE_DIM, typename OpBase>
1303 EntitiesFieldData::EntData &row_data) {
1305
1306 auto t_w = this->getFTensor0IntegrationWeight();
1307 auto t_base = row_data.getFTensor0N();
1308
1309 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1310 auto t_grad_y = getFTensor1FromMat<SPACE_DIM>(*yGradPtr);
1311
1313 const double alpha_constant = alphaConstant();
1314 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1315
1316 // get element volume
1317 const double vol = OpBase::getMeasure();
1318 const double c = (t_grad_y(i) * t_u(i)) * (t_w * vol * alpha_constant);
1319
1320 // get element volume
1321 int rr = 0;
1322 for (; rr != OpBase::nbRows; ++rr) {
1323 OpBase::locF[rr] += c * t_base;
1324 ++t_base;
1325 }
1326 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1327 ++t_base;
1328
1329 ++t_w; // move to another integration weight
1330 ++t_u;
1331 ++t_grad_y;
1332 }
1333
1335}
1336
1337template <int FIELD_DIM, int SPACE_DIM, typename OpBase>
1340 EntitiesFieldData::EntData &row_data) {
1342
1343 auto t_w = this->getFTensor0IntegrationWeight();
1344 auto t_base = row_data.getFTensor0N();
1345
1346 auto t_u = getFTensor1FromMat<SPACE_DIM>(*uPtr);
1347 auto t_grad_y = getFTensor2FromMat<FIELD_DIM, SPACE_DIM>(*yGradPtr);
1348
1351 const double alpha_constant = alphaConstant();
1352 for (size_t gg = 0; gg != OpBase::nbIntegrationPts; ++gg) {
1353
1354 // get element volume
1355 const double vol = OpBase::getMeasure();
1356
1358 t_c(J) = (t_grad_y(J, i) * t_u(i)) * (t_w * vol * alpha_constant);
1359
1360 auto t_nf = getFTensor1FromArray<FIELD_DIM, FIELD_DIM>(OpBase::locF);
1361 int rr = 0;
1362 for (; rr != OpBase::nbRows / FIELD_DIM; ++rr) {
1363 t_nf(J) += t_base * t_c(J);
1364 ++t_base;
1365 ++t_nf;
1366 }
1367 for (; rr < OpBase::nbRowBaseFunctions; ++rr)
1368 ++t_base;
1369 ++t_w; // move to another integration weight
1370
1371 ++t_u;
1372 ++t_grad_y;
1373 }
1374
1376}
1377
1378} // namespace MoFEM
1379
1380#endif // __LINEAR_FORMS_INTEGRATORS_HPP__
static Index< 'J', 3 > J
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 ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
CoordinateTypes
Coodinate system.
Definition: definitions.h:114
@ CYLINDRICAL
Definition: definitions.h:117
@ CARTESIAN
Definition: definitions.h:115
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
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
const double c
speed of light (cm/ns)
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: Common.hpp:10
boost::function< double(double)> TimeFun
Lambda function used to scale with time.
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 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)
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)
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)
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)