EdgeSlidingConstrains::CalculateEdgeBase Struct Reference

#include <users_modules/basic_finite_elements/src/SurfaceSlidingConstrains.hpp>

Static Public Member Functions
static MoFEMErrorCode	createTag (moab::Interface &moab, Tag &th0, Tag &th1, Tag &th2, Tag &th3)
static MoFEMErrorCode	numberSurfaces (moab::Interface &moab, Range edges, Range tris)
static MoFEMErrorCode	setTags (moab::Interface &moab, Range edges, Range tris, bool number_pathes=true, boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation > surface_orientation=nullptr, MoFEM::Interface *m_field_ptr=nullptr)
static MoFEMErrorCode	saveEdges (moab::Interface &moab, std::string name, Range edges, Range *faces=nullptr)

Detailed Description

Definition at line 641 of file SurfaceSlidingConstrains.hpp.

Member Function Documentation

createTag()

static MoFEMErrorCode EdgeSlidingConstrains::CalculateEdgeBase::createTag ( moab::Interface &  moab, Tag &  th0, Tag &  th1, Tag &  th2, Tag &  th3  )

inlinestatic

Definition at line 643 of file SurfaceSlidingConstrains.hpp.

                                                        {
      MoFEMFunctionBegin;
      double def_val[] = {0, 0, 0};
      CHKERR moab.tag_get_handle("EDGE_BASE0", 3, MB_TYPE_DOUBLE, th0,
                                 MB_TAG_CREAT | MB_TAG_SPARSE, def_val);
      CHKERR moab.tag_get_handle("EDGE_BASE1", 3, MB_TYPE_DOUBLE, th1,
                                 MB_TAG_CREAT | MB_TAG_SPARSE, def_val);
      int def_numb_val[] = {-1};
      CHKERR moab.tag_get_handle("PATCH_NUMBER", 1, MB_TYPE_INTEGER, th2,
                                 MB_TAG_CREAT | MB_TAG_SPARSE, def_numb_val);
      int def_orientation_val[] = {1};
      CHKERR moab.tag_get_handle("PATCH_ORIENTATION", 1, MB_TYPE_INTEGER, th3,
                                 MB_TAG_CREAT | MB_TAG_SPARSE,
                                 def_orientation_val);
 
      MoFEMFunctionReturn(0);
    }

numberSurfaces()

static MoFEMErrorCode EdgeSlidingConstrains::CalculateEdgeBase::numberSurfaces ( moab::Interface &  moab, Range  edges, Range  tris  )

inlinestatic

Definition at line 662 of file SurfaceSlidingConstrains.hpp.

                                                     {
      MoFEMFunctionBegin;
 
      auto get_edges = [&](const Range &ents) {
        Range edges;
        CHKERR moab.get_adjacencies(ents, 1, false, edges,
                                    moab::Interface::UNION);
        return edges;
      };
 
      auto get_face_adj = [&, get_edges](const Range &faces) {
        Range adj_faces;
        CHKERR moab.get_adjacencies(subtract(get_edges(faces), edges), 2, false,
                                    adj_faces, moab::Interface::UNION);
        return intersect(adj_faces, tris);
      };
 
      auto get_patch = [&, get_face_adj](const EntityHandle face) {
        Range patch_ents;
        patch_ents.insert(face);
        unsigned int nb0;
        do {
          nb0 = patch_ents.size();
          patch_ents.merge(get_face_adj(patch_ents));
        } while (nb0 != patch_ents.size());
        return patch_ents;
      };
 
      auto get_patches = [&]() {
        std::vector<Range> patches;
        while (!tris.empty()) {
          patches.push_back(get_patch(tris[0]));
          tris = subtract(tris, patches.back());
        }
        return patches;
      };
 
      Tag th0, th1, th2, th3;
      CHKERR createTag(moab, th0, th1, th2, th3);
 
      auto patches = get_patches();
      int pp = 0;
      for (auto patch : patches) {
        // cerr << "pp: " << pp << endl;
        // cerr << patch << endl;
        std::vector<int> tags_vals(patch.size(), pp);
        CHKERR moab.tag_set_data(th2, patch, &*tags_vals.begin());
        ++pp;
      }
 
      MoFEMFunctionReturn(0);
    }

saveEdges()

static MoFEMErrorCode EdgeSlidingConstrains::CalculateEdgeBase::saveEdges ( moab::Interface &  moab, std::string  name, Range  edges, Range *  faces = nullptr  )

inlinestatic

Definition at line 857 of file SurfaceSlidingConstrains.hpp.

                                                                         {
      MoFEMFunctionBegin;
      EntityHandle meshset;
      Tag ths[4];
      CHKERR createTag(moab, ths[0], ths[1], ths[2], ths[3]);
      CHKERR moab.create_meshset(MESHSET_SET, meshset);
      CHKERR moab.add_entities(meshset, edges);
      if (faces != nullptr) {
        CHKERR moab.add_entities(meshset, *faces);
      }
      CHKERR moab.write_file(name.c_str(), "VTK", "", &meshset, 1, ths, 4);
      CHKERR moab.delete_entities(&meshset, 1);
      MoFEMFunctionReturn(0);
    }

setTags()

static MoFEMErrorCode EdgeSlidingConstrains::CalculateEdgeBase::setTags ( moab::Interface &  moab, Range  edges, Range  tris, bool  number_pathes = true, boost::shared_ptr< SurfaceSlidingConstrains::DriverElementOrientation >  surface_orientation = nullptr, MoFEM::Interface *  m_field_ptr = nullptr  )

inlinestatic

Definition at line 716 of file SurfaceSlidingConstrains.hpp.

                                               {
      MoFEMFunctionBegin;
      Tag th0, th1, th2, th3;
      CHKERR createTag(moab, th0, th1, th2, th3);
      if (number_pathes) {
        CHKERR numberSurfaces(moab, edges, tris);
      }
      for (auto edge : edges) {
        Range adj_faces;
        CHKERR moab.get_adjacencies(&edge, 1, 2, false, adj_faces);
        adj_faces = intersect(adj_faces, tris);
        if (adj_faces.size() != 2) {
          SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
                   "Should be 2 faces adjacent to edge but is %d",
                   adj_faces.size());
        }
        VectorDouble3 v[2] = {VectorDouble3(3), VectorDouble3(3)};
        auto get_tensor_from_vec = [](VectorDouble3 &v) {
          return FTensor::Tensor1<FTensor::PackPtr<double *, 1>, 3>(
              &v[0], &v[1], &v[2]);
        };
 
        FTensor::Index<'i', 3> i;
        FTensor::Index<'j', 3> j;
        FTensor::Index<'k', 3> k;
 
        std::array<double, 3> areas;
        auto calculate_normals = [&, get_tensor_from_vec]() {
          int ff = 0;
          for (auto face : adj_faces) {
            double &x = (v[ff])[0];
            double &y = (v[ff])[1];
            double &z = (v[ff])[2];
            moab::Util::normal(&moab, face, x, y, z);
            auto t_n = get_tensor_from_vec(v[ff]);
            areas[ff] = sqrt(t_n(i) * t_n(i));
            t_n(i) /= areas[ff];
            int orientation;
            // FIXME: this tag is empty
            CHKERR moab.tag_get_data(th3, &face, 1, &orientation);
            if (orientation == -1) {
              t_n(i) *= -1;
            }
            if (surface_orientation && m_field_ptr) {
              CHKERR surface_orientation->getElementOrientation(*m_field_ptr,
                                                                face);
              int eo = surface_orientation->elementOrientation;
              if (eo == -1) {
                t_n(i) *= -1;
              }
            }
            ++ff;
          }
        };
        calculate_normals();
 
        auto get_patch_number = [&]() {
          std::vector<int> p = {0, 0};
          CHKERR moab.tag_get_data(th2, adj_faces, &*p.begin());
          return p;
        };
 
        std::array<EntityHandle, 2> faces_handles = {adj_faces[0],
                                                     adj_faces[1]};
        auto order_normals = [&, get_patch_number]() {
          auto p = get_patch_number();
          if (p[0] < p[1]) {
            v[0].swap(v[1]);
            faces_handles[0] = adj_faces[1];
            faces_handles[1] = adj_faces[0];
          }
        };
        order_normals();
 
        FTensor::Tensor1<double, 3> t_cross;
        auto t_n0 = get_tensor_from_vec(v[0]);
        auto t_n1 = get_tensor_from_vec(v[1]);
        t_cross(k) = FTensor::cross(t_n0(i), t_n1(j), k);
 
        auto get_base_for_coplanar_faces = [&]() {
          MoFEMFunctionBegin;
 
          std::vector<EntityHandle> face_conn;
          CHKERR moab.get_connectivity(&faces_handles[1], 1, face_conn, true);
          std::vector<EntityHandle> edge_conn;
          CHKERR moab.get_connectivity(&edge, 1, edge_conn, true);
          std::sort(face_conn.begin(), face_conn.end());
          std::sort(edge_conn.begin(), edge_conn.end());
          int n = 0;
          for (; n != 2; ++n)
            if (face_conn[n] != edge_conn[n])
              break;
          VectorDouble3 coords_face(3);
          CHKERR moab.get_coords(&face_conn[n], 1,
                                 &*coords_face.data().begin());
          VectorDouble6 coords_edge(6);
          CHKERR moab.get_coords(&edge_conn[0], 2,
                                 &*coords_edge.data().begin());
          auto t_edge_n0 = FTensor::Tensor1<double, 3>(
              coords_edge[0], coords_edge[1], coords_edge[2]);
          auto t_edge_n1 = FTensor::Tensor1<double, 3>(
              coords_edge[3], coords_edge[4], coords_edge[5]);
          auto t_face_n = get_tensor_from_vec(coords_face);
          t_face_n(i) -= t_edge_n0(i);
          FTensor::Tensor1<double, 3> t_delta;
          t_delta(i) = t_edge_n1(i) - t_edge_n0(i);
          t_delta(i) /= sqrt(t_delta(i) * t_delta(i));
          t_n1(k) = FTensor::cross(t_n0(i), t_delta(j), k);
          if (t_n1(i) * t_face_n(i) < 0)
            t_n1(i) *= -1;
 
          MoFEMFunctionReturn(0);
        };
 
        auto get_base_for_not_planar_faces = [&]() {
          MoFEMFunctionBeginHot;
          t_n1(k) = FTensor::cross(t_n0(i), t_cross(j), k);
          t_n1(i) /= sqrt(t_n1(i) * t_n1(i));
          MoFEMFunctionReturnHot(0);
        };
 
        // This a case when faces adjacent to edge are coplanar !!!
        constexpr double tol = 1e-6;
        if ((t_cross(i) * t_cross(i)) < tol)
          CHKERR get_base_for_coplanar_faces();
        else
          CHKERR get_base_for_not_planar_faces();
 
        VectorDouble3 &v0 = v[0];
        VectorDouble3 &v1 = v[1];
        CHKERR moab.tag_set_data(th0, &edge, 1, &v0[0]);
        CHKERR moab.tag_set_data(th1, &edge, 1, &v1[0]);
      }
      MoFEMFunctionReturn(0);
    }

---

The documentation for this struct was generated from the following file:

• users_modules/basic_finite_elements/src/SurfaceSlidingConstrains.hpp