50 #ifndef PANZER_FACE_TO_ELEMENT_IMPL_HPP
51 #define PANZER_FACE_TO_ELEMENT_IMPL_HPP
66 template <
typename LocalOrdinal,
typename GlobalOrdinal>
72 #ifndef PANZER_HIDE_DEPRECATED_CODE
77 template <
typename LocalOrdinal,
typename GlobalOrdinal>
85 template <
typename LocalOrdinal,
typename GlobalOrdinal>
90 initialize(conn, comm);
93 #ifndef PANZER_HIDE_DEPRECATED_CODE
98 template <
typename LocalOrdinal,
typename GlobalOrdinal>
104 initialize(conn, comm_world);
108 template <
typename LocalOrdinal,
typename GlobalOrdinal>
115 std::vector<std::string> block_ids;
121 std::vector<shards::CellTopology> ebt;
123 dimension = ebt[0].getDimension();
125 int my_rank = comm->getRank();
127 int nprocs = comm->getSize();
130 std::vector<GlobalOrdinal> element_GIDS;
131 for (
size_t iblk = 0 ; iblk < block_ids.size(); ++iblk) {
134 if ( dimension == 1 ) {
137 }
else if ( dimension == 2 ){
145 const std::vector<LocalOrdinal> &block_elems = conn.
getElementBlock(block_ids[iblk]);
146 for (
size_t i=0; i<block_elems.size(); ++i) {
148 element_GIDS.push_back(*connectivity);
157 std::vector<GlobalOrdinal> face_GIDS;
158 std::set<GlobalOrdinal> set_of_face_GIDS;
159 for (
size_t iblk = 0 ; iblk < block_ids.size(); ++iblk) {
162 if ( dimension == 1 ) {
165 }
else if ( dimension == 2 ){
173 const std::vector<LocalOrdinal> &block_elems = conn.
getElementBlock(block_ids[iblk]);
174 for (
size_t i=0; i<block_elems.size(); ++i) {
176 const panzer::GlobalOrdinal * connectivity = conn.
getConnectivity(block_elems[i]);
177 for (
int iface=0; iface<n_conn; ++iface)
178 set_of_face_GIDS.insert(connectivity[iface]);
181 face_GIDS.insert(face_GIDS.begin(), set_of_face_GIDS.begin(), set_of_face_GIDS.end());
184 owned_face_map = Tpetra::createOneToOne(face_map);
193 face2elem_mv->putScalar(-1-shift);
195 auto f2e = face2elem_mv->getLocalViewHost(Tpetra::Access::ReadWrite);
196 auto b1 = Kokkos::subview(f2e,Kokkos::ALL(),0);
197 auto e1 = Kokkos::subview(f2e,Kokkos::ALL(),1);
198 auto p1 = Kokkos::subview(f2e,Kokkos::ALL(),2);
199 auto l1 = Kokkos::subview(f2e,Kokkos::ALL(),3);
200 auto b2 = Kokkos::subview(f2e,Kokkos::ALL(),4);
201 auto e2 = Kokkos::subview(f2e,Kokkos::ALL(),5);
202 auto p2 = Kokkos::subview(f2e,Kokkos::ALL(),6);
203 auto l2 = Kokkos::subview(f2e,Kokkos::ALL(),7);
205 GlobalOrdinal my_elem = 0;
206 for (
size_t iblk = 0 ; iblk < block_ids.size(); ++iblk) {
209 if ( dimension == 1 ) {
212 }
else if ( dimension == 2 ){
220 const std::vector<LocalOrdinal> &block_elems = conn.
getElementBlock(block_ids[iblk]);
221 for (
size_t i=0; i<block_elems.size(); ++i) {
223 const panzer::GlobalOrdinal * connectivity = conn.
getConnectivity(block_elems[i]);
224 for (
int iface=0; iface<n_conn; ++iface) {
225 LocalOrdinal f = face_map->getLocalElement(connectivity[iface]);
231 e1[f] = elem_map->getGlobalElement(my_elem)+shift;
232 p1[f] = my_rank+shift;
234 }
else if (b2[f] < 0){
236 e2[f] = elem_map->getGlobalElement(my_elem)+shift;
237 p2[f] = my_rank+shift;
248 Import imp(owned_face_map, face_map);
249 Export exp(face_map, owned_face_map);
250 owned_face2elem_mv->doExport(*face2elem_mv, exp, Tpetra::ADD);
253 auto f2e = face2elem_mv->getLocalViewHost(Tpetra::Access::ReadWrite);
254 auto b1 = Kokkos::subview(f2e,Kokkos::ALL(),0);
255 auto e1 = Kokkos::subview(f2e,Kokkos::ALL(),1);
256 auto p1 = Kokkos::subview(f2e,Kokkos::ALL(),2);
257 auto l1 = Kokkos::subview(f2e,Kokkos::ALL(),3);
258 auto b2 = Kokkos::subview(f2e,Kokkos::ALL(),4);
259 auto e2 = Kokkos::subview(f2e,Kokkos::ALL(),5);
260 auto p2 = Kokkos::subview(f2e,Kokkos::ALL(),6);
261 auto l2 = Kokkos::subview(f2e,Kokkos::ALL(),7);
263 auto of2e = owned_face2elem_mv->getLocalViewHost(Tpetra::Access::ReadWrite);
264 auto ob1 = Kokkos::subview(of2e,Kokkos::ALL(),0);
265 auto oe1 = Kokkos::subview(of2e,Kokkos::ALL(),1);
266 auto op1 = Kokkos::subview(of2e,Kokkos::ALL(),2);
267 auto ol1 = Kokkos::subview(of2e,Kokkos::ALL(),3);
268 auto ob2 = Kokkos::subview(of2e,Kokkos::ALL(),4);
269 auto oe2 = Kokkos::subview(of2e,Kokkos::ALL(),5);
270 auto op2 = Kokkos::subview(of2e,Kokkos::ALL(),6);
271 auto ol2 = Kokkos::subview(of2e,Kokkos::ALL(),7);
276 for (
size_t i=0; i<ob1.size();++i){
279 assert(b1[i] >= shift);
281 LocalOrdinal shared_local_id = face_map->getLocalElement(owned_face_map->getGlobalElement(i));
283 if (ob1[i] == b1[shared_local_id] && ob2[i] == b2[shared_local_id] &&
284 oe1[i] == e1[shared_local_id] && oe2[i] == e2[shared_local_id]) {
291 if (ob1[i] < b1[shared_local_id] || oe1[i] < e1[shared_local_id]) {
295 if ( ob1[i] > b1[shared_local_id] || oe1[i] > e1[shared_local_id]) {
297 assert(ob2[i] < 0 && oe2[i] < 0);
298 ob2[i] = ob1[i] - b1[shared_local_id];
299 oe2[i] = oe1[i] - e1[shared_local_id];
300 op2[i] = op1[i] - p1[shared_local_id];
301 ol2[i] = ol1[i] - l1[shared_local_id];
303 ob1[i] = b1[shared_local_id];
304 oe1[i] = e1[shared_local_id];
305 op1[i] = p1[shared_local_id];
306 ol1[i] = l1[shared_local_id];
308 assert(op1[i] >=0 && op2[i] >= 0 && op1[i] < nprocs+shift && op2[i] < nprocs+shift);
312 face2elem_mv->doImport(*owned_face2elem_mv, imp, Tpetra::REPLACE);
317 face_map_ = face_map;
318 LocalOrdinal nfaces = face_map_->getLocalNumElements();
319 elems_by_face_ = PHX::View<GlobalOrdinal *[2]>(
"FaceToElement::elems_by_face_", nfaces);
320 lidx_by_face_ = PHX::View<int *[2]>(
"FaceToElement::elems_by_face_", nfaces);
321 blocks_by_face_ = PHX::View<int *[2]> (
"FaceToElement::blocks_by_face_", nfaces);
322 procs_by_face_ = PHX::View<int *[2]> (
"FaceToElement::procs_by_face_", nfaces);
323 auto elems_by_face_h = Kokkos::create_mirror_view(elems_by_face_);
324 auto blocks_by_face_h = Kokkos::create_mirror_view(blocks_by_face_);
325 auto procs_by_face_h = Kokkos::create_mirror_view(procs_by_face_);
326 auto lidx_by_face_h = Kokkos::create_mirror_view(lidx_by_face_);
328 auto f2e = face2elem_mv->getLocalViewHost(Tpetra::Access::ReadWrite);
329 auto b1 = Kokkos::subview(f2e,Kokkos::ALL(),0);
330 auto e1 = Kokkos::subview(f2e,Kokkos::ALL(),1);
331 auto p1 = Kokkos::subview(f2e,Kokkos::ALL(),2);
332 auto l1 = Kokkos::subview(f2e,Kokkos::ALL(),3);
333 auto b2 = Kokkos::subview(f2e,Kokkos::ALL(),4);
334 auto e2 = Kokkos::subview(f2e,Kokkos::ALL(),5);
335 auto p2 = Kokkos::subview(f2e,Kokkos::ALL(),6);
336 auto l2 = Kokkos::subview(f2e,Kokkos::ALL(),7);
338 for (LocalOrdinal i=0; i< nfaces; ++i) {
339 elems_by_face_h (i,0) = e1[i]-shift;
340 elems_by_face_h (i,1) = e2[i]-shift;
341 blocks_by_face_h(i,0) = b1[i]-shift;
342 blocks_by_face_h(i,1) = b2[i]-shift;
343 procs_by_face_h (i,0) = p1[i]-shift;
344 procs_by_face_h (i,1) = p2[i]-shift;
345 lidx_by_face_h (i,0) = l1[i]-shift;
346 lidx_by_face_h (i,1) = l2[i]-shift;
348 if(elems_by_face_h (i,0) < 0){
349 elems_by_face_h (i,0) = -1;
351 if(elems_by_face_h (i,1) < 0){
352 elems_by_face_h (i,1) = -1;
355 Kokkos::deep_copy(elems_by_face_, elems_by_face_h);
356 Kokkos::deep_copy(blocks_by_face_, blocks_by_face_h);
357 Kokkos::deep_copy(procs_by_face_, procs_by_face_h);
358 Kokkos::deep_copy(lidx_by_face_, lidx_by_face_h);
Tpetra::Import< LocalOrdinal, GlobalOrdinal, NodeType > Import
virtual LocalOrdinal getConnectivitySize(LocalOrdinal localElmtId) const =0
void initialize(panzer::ConnManager &conn)
Tpetra::Map< LocalOrdinal, GlobalOrdinal, NodeType > Map
virtual const std::vector< LocalOrdinal > & getElementBlock(const std::string &blockID) const =0
virtual void getElementBlockIds(std::vector< std::string > &elementBlockIds) const =0
Pure virtual base class for supplying mesh connectivity information to the DOF Manager.
virtual void getElementBlockTopologies(std::vector< shards::CellTopology > &elementBlockTopologies) const =0
Tpetra::Export< LocalOrdinal, GlobalOrdinal, NodeType > Export
Tpetra::MultiVector< GlobalOrdinal, LocalOrdinal, GlobalOrdinal, NodeType > GOMultiVector
virtual void buildConnectivity(const FieldPattern &fp)=0
virtual const GlobalOrdinal * getConnectivity(LocalOrdinal localElmtId) const =0