1 #ifndef FILE_IMPROVE2
2 #define FILE_IMPROVE2
3
AppendEdges(const Element2d & elem,PointIndex pi,Array<std::tuple<PointIndex,PointIndex>> & edges)4 inline void AppendEdges( const Element2d & elem, PointIndex pi, Array<std::tuple<PointIndex,PointIndex>> & edges )
5 {
6 for (int j = 0; j < 3; j++)
7 {
8 PointIndex pi0 = elem[j];
9 PointIndex pi1 = elem[(j+1)%3];
10 if (pi1 < pi0) Swap(pi0, pi1);
11 if(pi0==pi)
12 edges.Append(std::make_tuple(pi0, pi1));
13 }
14 }
15
AppendEdges(const Element & elem,PointIndex pi,Array<std::tuple<PointIndex,PointIndex>> & edges)16 inline void AppendEdges( const Element & elem, PointIndex pi, Array<std::tuple<PointIndex,PointIndex>> & edges )
17 {
18 static constexpr int tetedges[6][2] =
19 { { 0, 1 }, { 0, 2 }, { 0, 3 },
20 { 1, 2 }, { 1, 3 }, { 2, 3 } };
21
22 if(elem.flags.fixed)
23 return;
24 for (int j = 0; j < 6; j++)
25 {
26 PointIndex pi0 = elem[tetedges[j][0]];
27 PointIndex pi1 = elem[tetedges[j][1]];
28 if (pi1 < pi0) Swap(pi0, pi1);
29 if(pi0==pi)
30 edges.Append(std::make_tuple(pi0, pi1));
31 }
32 }
33
34 template<typename TINDEX>
BuildEdgeList(const Mesh & mesh,const Table<TINDEX,PointIndex> & elementsonnode,Array<std::tuple<PointIndex,PointIndex>> & edges)35 void BuildEdgeList( const Mesh & mesh, const Table<TINDEX, PointIndex> & elementsonnode, Array<std::tuple<PointIndex, PointIndex>> & edges )
36 {
37 static_assert(is_same_v<TINDEX, ElementIndex>||is_same_v<TINDEX,SurfaceElementIndex>, "Invalid type for TINDEX");
38 static Timer tbuild_edges("Build edges"); RegionTimer reg(tbuild_edges);
39
40 int ntasks = 4*ngcore::TaskManager::GetMaxThreads();
41 Array<Array<std::tuple<PointIndex,PointIndex>>> task_edges(ntasks);
42
43 ParallelFor(IntRange(ntasks), [&] (int ti)
44 {
45 auto myrange = mesh.Points().Range().Split(ti, ntasks);
46 ArrayMem<std::tuple<PointIndex,PointIndex>, 100> local_edges;
47 for (auto pi : myrange)
48 {
49 local_edges.SetSize(0);
50
51 for(auto ei : elementsonnode[pi])
52 {
53 const auto & elem = mesh[ei];
54 if (elem.IsDeleted()) continue;
55
56 AppendEdges(elem, pi, local_edges);
57 }
58 QuickSort(local_edges);
59
60 auto edge_prev = std::make_tuple<PointIndex, PointIndex>(-1,-1);
61
62 for(auto edge : local_edges)
63 if(edge != edge_prev)
64 {
65 task_edges[ti].Append(edge);
66 edge_prev = edge;
67 }
68 }
69 }, ntasks);
70
71 int num_edges = 0;
72 for (auto & edg : task_edges)
73 num_edges += edg.Size();
74 edges.SetAllocSize(num_edges);
75 for (auto & edg : task_edges)
76 edges.Append(edg);
77 }
78
79
80 class Neighbour
81 {
82 int nr[3];
83 int orient[3];
84
85 public:
Neighbour()86 Neighbour () { ; }
87
SetNr(int side,int anr)88 void SetNr (int side, int anr) { nr[side] = anr; }
GetNr(int side)89 int GetNr (int side) { return nr[side]; }
90
SetOrientation(int side,int aorient)91 void SetOrientation (int side, int aorient) { orient[side] = aorient; }
GetOrientation(int side)92 int GetOrientation (int side) { return orient[side]; }
93 };
94
95 ///
96 class MeshOptimize2d
97 {
98 int faceindex = 0;
99 int improveedges = 0;
100 double metricweight = 0.;
101 int writestatus = 1;
102 Mesh& mesh;
103 const NetgenGeometry& geo;
104 public:
105 ///
MeshOptimize2d(Mesh & amesh)106 MeshOptimize2d(Mesh& amesh) : mesh(amesh), geo(*mesh.GetGeometry())
107 {}
~MeshOptimize2d()108 virtual ~MeshOptimize2d() { ; }
109 ///
110 void ImproveMesh (const MeshingParameters & mp);
111 void ImproveMeshJacobian (const MeshingParameters & mp);
112 void ImproveVolumeMesh ();
113 void ProjectBoundaryPoints(NgArray<int> & surfaceindex,
114 const NgArray<Point<3>* > & from, NgArray<Point<3>* > & dest);
115
116 bool EdgeSwapping (const int usemetric, Array<Neighbour> &neighbors, Array<bool> &swapped,
117 const SurfaceElementIndex t1, const int edge, const int t, Array<int,PointIndex> &pdef, const bool check_only=false);
118 void EdgeSwapping (int usemetric);
119 void CombineImprove ();
120 void SplitImprove ();
121
122 void GenericImprove ();
123
124
SetFaceIndex(int fi)125 void SetFaceIndex (int fi) { faceindex = fi; }
SetImproveEdges(int ie)126 void SetImproveEdges (int ie) { improveedges = ie; }
SetMetricWeight(double mw)127 void SetMetricWeight (double mw) { metricweight = mw; }
SetWriteStatus(int ws)128 void SetWriteStatus (int ws) { writestatus = ws; }
129
130
131 /// liefert zu einem 3d-Punkt die geominfo (Dreieck) und liefert 1, wenn erfolgreich,
132 /// 0, wenn nicht (Punkt ausserhalb von chart)
133 ///
134
135 void CheckMeshApproximation (Mesh & mesh);
136
137
138 ///
139 friend class Opti2SurfaceMinFunction;
140 ///
141 friend class Opti2EdgeMinFunction;
142 ///
143 friend double Opti2FunctionValueGrad (const Vector & x, Vector & grad);
144 ///
145 friend double Opti2EdgeFunctionValueGrad (const Vector & x, Vector & grad);
146
147
148
149 };
150
151
152 extern void CalcTriangleBadness (double x2, double x3, double y3,
153 double metricweight,
154 double h, double & badness,
155 double & g1x, double & g1y);
156
157
158
159
160 extern double CalcTriangleBadness (const Point<3> & p1,
161 const Point<3> & p2,
162 const Point<3> & p3,
163 double metricweight,
164 double h);
165
166 extern double CalcTriangleBadness (const Point<3> & p1,
167 const Point<3> & p2,
168 const Point<3> & p3,
169 const Vec<3> & n,
170 double metricweight,
171 double h);
172
173 #endif
174
175
176