/[escript]/trunk/finley/src/Mesh_joinFaces.cpp
ViewVC logotype

Contents of /trunk/finley/src/Mesh_joinFaces.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6939 - (show annotations)
Mon Jan 20 03:37:18 2020 UTC (4 months, 1 week ago) by uqaeller
File size: 5218 byte(s)
Updated the copyright header.


1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2020 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014-2017 by Centre for Geoscience Computing (GeoComp)
14 * Development from 2019 by School of Earth and Environmental Sciences
15 **
16 *****************************************************************************/
17
18
19 /****************************************************************************
20
21 Finley: Mesh
22
23 detects faces in the mesh that match and replaces it by step elements
24
25 *****************************************************************************/
26
27 #include "FinleyDomain.h"
28
29 #include <escript/index.h>
30
31 namespace finley {
32
33 void FinleyDomain::joinFaces(double safety_factor, double tolerance, bool optimize)
34 {
35 if (m_mpiInfo->size > 1) {
36 throw escript::NotImplementedError("Mesh::joinFaces: MPI is not supported yet.");
37 }
38 if (!m_contactElements) {
39 throw escript::ValueError("Mesh::joinFaces: No contact elements present.");
40 }
41 if (!m_faceElements)
42 return;
43
44 const_ReferenceElement_ptr faceRefElement(m_faceElements->referenceElementSet->borrowReferenceElement(false));
45 const_ReferenceElement_ptr contactRefElement(m_contactElements->referenceElementSet->borrowReferenceElement(false));
46
47 if (faceRefElement->Type->numNodesOnFace <= 0) {
48 std::stringstream ss;
49 ss << "Mesh::joinFaces: joining faces cannot be applied to face "
50 "elements of type " << faceRefElement->Type->Name;
51 throw escript::ValueError(ss.str());
52 }
53
54 if (contactRefElement->Type->numNodes != 2*faceRefElement->Type->numNodes) {
55 std::stringstream ss;
56 ss << "Mesh::joinFaces: contact element file for "
57 << contactRefElement->Type->Name << " needs to hold elements "
58 "created from face elements " << faceRefElement->Type->Name;
59 throw escript::ValueError(ss.str());
60 }
61
62 const int NN = m_faceElements->numNodes;
63 const int NN_Contact = m_contactElements->numNodes;
64
65 // allocate work arrays
66 int* elem1 = new int[m_faceElements->numElements];
67 int* elem0 = new int[m_faceElements->numElements];
68 index_t* elem_mask = new index_t[m_faceElements->numElements];
69 int* matching_nodes_in_elem1 = new int[m_faceElements->numElements*NN];
70
71 // find the matching face elements
72 int numPairs;
73 findMatchingFaces(safety_factor, tolerance, &numPairs, elem0, elem1, matching_nodes_in_elem1);
74 // get a list of the face elements to be kept
75 #pragma omp parallel for
76 for (index_t e = 0; e < m_faceElements->numElements; e++)
77 elem_mask[e] = 1;
78
79 for (int e = 0; e < numPairs; e++) {
80 elem_mask[elem0[e]] = 0;
81 elem_mask[elem1[e]] = 0;
82 }
83 dim_t new_numFaceElements = 0;
84 // OMP
85 for (index_t e = 0; e < m_faceElements->numElements; e++) {
86 if (elem_mask[e] > 0) {
87 elem_mask[new_numFaceElements] = e;
88 new_numFaceElements++;
89 }
90 }
91 // allocate new face element and Contact element files
92 ElementFile* newContactElementsFile = new ElementFile(m_contactElements->referenceElementSet, m_mpiInfo);
93 ElementFile* newFaceElementsFile = new ElementFile(m_faceElements->referenceElementSet, m_mpiInfo);
94 newContactElementsFile->allocTable(numPairs+m_contactElements->numElements);
95 newFaceElementsFile->allocTable(new_numFaceElements);
96 // copy the old elements over
97 // get the face elements which are still in use
98 newFaceElementsFile->gather(elem_mask, m_faceElements);
99 // get the contact elements which are still in use
100 newContactElementsFile->copyTable(0, 0, 0, m_contactElements);
101 dim_t c = m_contactElements->numElements;
102 // OMP
103 for (int e = 0; e < numPairs; e++) {
104 const int e0 = elem0[e];
105 const int e1 = elem1[e];
106 newContactElementsFile->Id[c] = std::min(m_faceElements->Id[e0], m_faceElements->Id[e1]);
107 newContactElementsFile->Tag[c] = std::min(m_faceElements->Tag[e0], m_faceElements->Tag[e1]);
108 newContactElementsFile->Color[c] = e;
109 for (int i = 0; i < NN; i++)
110 newContactElementsFile->Nodes[INDEX2(i,c,NN_Contact)] =
111 m_faceElements->Nodes[INDEX2(i,e0,NN)];
112 for (int i = 0; i < NN; i++)
113 newContactElementsFile->Nodes[INDEX2(i+NN,c,NN_Contact)] =
114 matching_nodes_in_elem1[INDEX2(i,e,NN)];
115 c++;
116 }
117 newContactElementsFile->minColor = 0;
118 newContactElementsFile->maxColor = numPairs-1;
119 // set new face and Contact elements
120
121 delete m_faceElements;
122 m_faceElements = newFaceElementsFile;
123
124 delete m_contactElements;
125 m_contactElements = newContactElementsFile;
126
127 prepare(optimize);
128
129 delete[] elem1;
130 delete[] elem0;
131 delete[] matching_nodes_in_elem1;
132 delete[] elem_mask;
133 }
134
135 } // namespace finley
136

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision
svn:mergeinfo /branches/4.0fordebian/finley/src/Mesh_joinFaces.cpp:5567-5588 /branches/lapack2681/finley/src/Mesh_joinFaces.cpp:2682-2741 /branches/pasowrap/finley/src/Mesh_joinFaces.cpp:3661-3674 /branches/py3_attempt2/finley/src/Mesh_joinFaces.cpp:3871-3891 /branches/restext/finley/src/Mesh_joinFaces.cpp:2610-2624 /branches/ripleygmg_from_3668/finley/src/Mesh_joinFaces.cpp:3669-3791 /branches/stage3.0/finley/src/Mesh_joinFaces.cpp:2569-2590 /branches/symbolic_from_3470/finley/src/Mesh_joinFaces.cpp:3471-3974 /branches/trilinos_from_5897/finley/src/Mesh_joinFaces.cpp:5898-6118 /release/3.0/finley/src/Mesh_joinFaces.cpp:2591-2601 /release/4.0/finley/src/Mesh_joinFaces.cpp:5380-5406 /trunk/finley/src/Mesh_joinFaces.cpp:4257-4344

  ViewVC Help
Powered by ViewVC 1.1.26