v0.15.0
Loading...
Searching...
No Matches
SeepageOps.hpp
Go to the documentation of this file.
1/* This file is part of MoFEM.
2 * MoFEM is free software: you can redistribute it and/or modify it under
3 * the terms of the GNU Lesser General Public License as published by the
4 * Free Software Foundation, either version 3 of the License, or (at your
5 * option) any later version.
6 *
7 * MoFEM is distributed in the hope that it will be useful, but WITHOUT
8 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
10 * License for more details.
11 *
12 * You should have received a copy of the GNU Lesser General Public
13 * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
14
15/** \file SeepageOps.hpp
16 * \example SeepageOps.hpp
17 */
18
19namespace SeepageOps {
20
21template <int DIM>
23 : public AssemblyDomainEleOp { // changed opfaceele to AssemblyDomainEleOp
24public:
25 OpDomainRhsHydrostaticStress(std::string field_name1,
26 boost::shared_ptr<VectorDouble> h_ptr,
27 double specific_weight_water = 9.81)
28 : AssemblyDomainEleOp(field_name1, field_name1, DomainEleOp::OPROW),
29 hPtr(h_ptr), specificWeightWater(specific_weight_water) {}
30
31 MoFEMErrorCode iNtegrate(DataForcesAndSourcesCore::EntData &data) {
33
34 const int nb_dofs = data.getIndices().size();
35
36 if (nb_dofs) {
37 // locRhs.resize(nb_dofs, false);
38 // locRhs.clear();
39 auto &nf = AssemblyDomainEleOp::locF;
40 // get element area
41 const double area = getMeasure();
42
43 // get number of integration points
44 const int nb_integration_points = getGaussPts().size2();
45 // get integration weights
46 auto t_w = getFTensor0IntegrationWeight();
47
48 // get base function
49 auto t_base_diff = data.getFTensor1DiffN<DIM>();
50
51 FTensor::Index<'i', DIM> i;
52
53 auto t_h = getFTensor0FromVec(*hPtr);
54 for (int gg = 0; gg != nb_integration_points; gg++) {
55 auto t_nf = getFTensor1FromPtr<DIM>(&nf[0]);
56
57 const double a = t_w * area * specificWeightWater * t_h;
58
59 for (int rr = 0; rr != nb_dofs / DIM; rr++) {
60 t_nf(i) -= t_base_diff(i) * a;
61
62 // move to the next base function
63 ++t_base_diff; // moves the pointer to the next shape function
64 ++t_nf;
65 }
66
67 // move to the weight of the next integration point
68 ++t_w;
69 ++t_h;
70 }
71 }
72
74 }
75
76private:
77 // VectorDouble locRhs;
79 boost::shared_ptr<VectorDouble> hPtr;
80};
81
82} // namespace SeepageOps
constexpr double a
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
FTensor::Index< 'i', SPACE_DIM > i
MoFEMErrorCode iNtegrate(DataForcesAndSourcesCore::EntData &data)
OpDomainRhsHydrostaticStress(std::string field_name1, boost::shared_ptr< VectorDouble > h_ptr, double specific_weight_water=9.81)
boost::shared_ptr< VectorDouble > hPtr