MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
transfermap.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #include "submesh.hpp"
13 #include "transfermap.hpp"
14 #include "submesh_utils.hpp"
15 
16 using namespace mfem;
17 
19  const GridFunction &dst)
20 {
21  const FiniteElementSpace *parentfes = nullptr, *subfes1 = nullptr,
22  *subfes2 = nullptr;
23 
24  if (SubMesh::IsSubMesh(src.FESpace()->GetMesh()) &&
26  {
27  SubMesh* src_sm = static_cast<SubMesh*>(src.FESpace()->GetMesh());
28  SubMesh* dst_sm = static_cast<SubMesh*>(dst.FESpace()->GetMesh());
29 
30  // There is no immediate relation and both src and dst come from a
31  // SubMesh, check if they have an equivalent root parent.
32  if (SubMeshUtils::GetRootParent<SubMesh, Mesh>(*src_sm) !=
33  SubMeshUtils::GetRootParent<SubMesh, Mesh>(*dst_sm))
34  {
35  MFEM_ABORT("Can't find a relation between the two GridFunctions");
36  }
37 
39 
40  root_fes_ = new FiniteElementSpace(*src.FESpace(),
41  const_cast<Mesh *>(SubMeshUtils::GetRootParent<SubMesh, Mesh>
42  (*src_sm)));
43  subfes1 = src.FESpace();
44  subfes2 = dst.FESpace();
45 
47  *root_fes_,
48  src_sm->GetFrom(),
49  src_sm->GetParentElementIDMap(),
50  sub1_to_parent_map_);
51 
53  *root_fes_,
54  dst_sm->GetFrom(),
55  dst_sm->GetParentElementIDMap(),
56  sub2_to_parent_map_);
57 
58  z_.SetSize(root_fes_->GetVSize());
59  }
60  else if (SubMesh::IsSubMesh(src.FESpace()->GetMesh()))
61  {
63  SubMesh* src_sm = static_cast<SubMesh*>(src.FESpace()->GetMesh());
64  subfes1 = src.FESpace();
65  parentfes = dst.FESpace();
67  *parentfes,
68  src_sm->GetFrom(),
69  src_sm->GetParentElementIDMap(),
70  sub1_to_parent_map_);
71  }
72  else if (SubMesh::IsSubMesh(dst.FESpace()->GetMesh()))
73  {
75  SubMesh* dst_sm = static_cast<SubMesh*>(dst.FESpace()->GetMesh());
76  subfes1 = dst.FESpace();
77  parentfes = src.FESpace();
79  *parentfes,
80  dst_sm->GetFrom(),
81  dst_sm->GetParentElementIDMap(),
82  sub1_to_parent_map_);
83  }
84  else
85  {
86  MFEM_ABORT("Trying to do a transfer between GridFunctions but none of them is defined on a SubMesh");
87  }
88 }
89 
91  GridFunction &dst) const
92 {
93  if (category_ == TransferCategory::ParentToSubMesh)
94  {
95  // dst = S1^T src
96  for (int i = 0; i < sub1_to_parent_map_.Size(); i++)
97  {
98  dst(i) = src(sub1_to_parent_map_[i]);
99  }
100  }
101  else if (category_ == TransferCategory::SubMeshToParent)
102  {
103  // dst = G S1 src
104  // = G z
105  //
106  // G is identity if the partitioning matches
107 
108  for (int i = 0; i < sub1_to_parent_map_.Size(); i++)
109  {
110  dst(sub1_to_parent_map_[i]) = src(i);
111  }
112  }
113  else if (category_ == TransferCategory::SubMeshToSubMesh)
114  {
115  // dst = S2^T G (S1 src (*) S2 dst)
116  //
117  // G is identity if the partitioning matches
118 
119  z_ = 0.0;
120 
121  for (int i = 0; i < sub2_to_parent_map_.Size(); i++)
122  {
123  z_(sub2_to_parent_map_[i]) = dst(i);
124  }
125 
126  for (int i = 0; i < sub1_to_parent_map_.Size(); i++)
127  {
128  z_(sub1_to_parent_map_[i]) = src(i);
129  }
130 
131  for (int i = 0; i < sub2_to_parent_map_.Size(); i++)
132  {
133  dst(i) = z_(sub2_to_parent_map_[i]);
134  }
135  }
136  else
137  {
138  MFEM_ABORT("unknown TransferCategory: " << category_);
139  }
140 }
141 
143 {
144  delete root_fes_;
145 }
int Size() const
Return the logical size of the array.
Definition: array.hpp:138
int GetVSize() const
Return the number of vector dofs, i.e. GetNDofs() x GetVDim().
Definition: fespace.hpp:587
void BuildVdofToVdofMap(const FiniteElementSpace &subfes, const FiniteElementSpace &parentfes, const SubMesh::From &from, const Array< int > &parent_element_ids, Array< int > &vdof_to_vdof_map)
Build the vdof to vdof mapping between two FiniteElementSpace objects.
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
const Array< int > & GetParentElementIDMap() const
Get the parent element id map.
Definition: submesh.hpp:108
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:513
void Transfer(const GridFunction &src, GridFunction &dst) const
Transfer the source GridFunction to the destination GridFunction.
Definition: transfermap.cpp:90
Mesh * GetMesh() const
Returns the mesh.
Definition: fespace.hpp:441
static bool IsSubMesh(const Mesh *m)
Check if Mesh m is a SubMesh.
Definition: submesh.hpp:162
FiniteElementSpace * FESpace()
Definition: gridfunc.hpp:668
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:96
From GetFrom() const
Get the From indicator.
Definition: submesh.hpp:98
Subdomain representation of a topological parent in another Mesh.
Definition: submesh.hpp:42
TransferMap(const GridFunction &src, const GridFunction &dst)
Construct a new TransferMap object which transfers degrees of freedom from the source GridFunction to...
Definition: transfermap.cpp:18