From bb99e596d34a9256c4d1d2d235eb804b39c033ab Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?=D0=A1=D0=B5=D1=80=D0=B3=D0=B5=D0=B9=20=D0=9D=D0=B8=D0=BA?=
 =?UTF-8?q?=D1=83=D0=BB=D0=BA=D0=B8=D0=BD?= <sa.nikulkin@gmail.com>
Date: Wed, 1 Dec 2021 12:43:01 +0000
Subject: [PATCH] Update MeshLoader.cpp

---
 src/MeshLoader.cpp | 69 ++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 69 insertions(+)

diff --git a/src/MeshLoader.cpp b/src/MeshLoader.cpp
index fc4c200..ec0e39b 100644
--- a/src/MeshLoader.cpp
+++ b/src/MeshLoader.cpp
@@ -4,6 +4,8 @@
 #include <list>
 #include <set>
 #include <unordered_set>
+#include <ranges>
+#include <cmath>
 
 #include "MeshLoader.h"
 
@@ -220,3 +222,70 @@ std::vector<std::vector<Node>> MeshLoader::NodesNeighbors() {
 
     return ans;
 }
+
+double GetAngleWorker(Node& a, Node& b, Node& c, Node& d) {
+    double tmp1x = (c.coord[1] - a.coord[1]) * (b.coord[2] - a.coord[2]) - 
+                  (b.coord[1] - a.coord[1]) * (c.coord[2] - a.coord[2]);
+
+    double tmp1y = (b.coord[0] - a.coord[0]) * (c.coord[2] - a.coord[2]) -
+                  (c.coord[0] - a.coord[0]) * (b.coord[2] - a.coord[2]);
+
+    double tmp1z = (c.coord[0] - a.coord[0]) * (b.coord[1] - a.coord[1]) -
+                  (b.coord[0] - a.coord[0]) * (c.coord[1] - a.coord[1]);
+
+    double tmp1 = tmp1x * tmp1x + tmp1y * tmp1y + tmp1z * tmp1z;
+
+    double tmp2x = (d.coord[1] - a.coord[1]) * (b.coord[2] - a.coord[2]) -
+                   (b.coord[1] - a.coord[1]) * (d.coord[2] - a.coord[2]);
+
+    double tmp2y = (b.coord[0] - a.coord[0]) * (d.coord[2] - a.coord[2]) -
+                   (d.coord[0] - a.coord[0]) * (b.coord[2] - a.coord[2]);
+
+    double tmp2z = (d.coord[0] - a.coord[0]) * (b.coord[1] - a.coord[1]) -
+                   (b.coord[0] - a.coord[0]) * (d.coord[1] - a.coord[1]);
+
+    double tmp2 = tmp2x * tmp2x + tmp2y * tmp2y + tmp2z * tmp2z;
+
+    double chsl = tmp1 + tmp2 - ((d.coord[0] - c.coord[0]) * (d.coord[0] - c.coord[0]) +
+                                 (d.coord[1] - c.coord[1]) * (d.coord[1] - c.coord[1]) +
+                                 (d.coord[2] - c.coord[2]) * (d.coord[2] - c.coord[2])) * 
+                               ((b.coord[0] - a.coord[0]) * (b.coord[0] - a.coord[0]) +
+                                (b.coord[1] - a.coord[1]) * (b.coord[1] - a.coord[1]) +
+                                (b.coord[2] - a.coord[2]) * (b.coord[2] - a.coord[2]));
+    
+    double znm = 2 * std::sqrt(tmp1) * std::sqrt(tmp2);
+    
+    return acos(chsl / znm);
+}
+
+double MeshLoader::GetAngle(const FiniteElement& elem) {
+    if (elem.nodes_id.size() < 4)
+        return 0.0;
+    Node a = GetNodeById(elem.nodes_id[0]);
+    Node b = GetNodeById(elem.nodes_id[1]);
+    Node c = GetNodeById(elem.nodes_id[2]);
+    Node d = GetNodeById(elem.nodes_id[3]);
+
+    double alpha = GetAngleWorker(a, b, c, d);
+    double betha = GetAngleWorker(b, c, a, d);
+    double hama = GetAngleWorker(a, c, b, d);
+    
+    return std::max(std::max(alpha, betha), hama);
+}
+
+void MeshLoader::SortByAngle(const std::string& filepath) {
+    std::vector<FiniteElement> elements = FiniteElements();
+    if (elements.empty())
+        return;
+    std::ofstream file(filepath + ".txt");
+    std::ranges::sort(elements, std::less<>{}, [&](const FiniteElement& elem) {
+        return GetAngle(elem);
+        });
+    std::ranges::copy(
+        elements | std::ranges::views::drop(9) | std::ranges::views::take(10) |
+        std::ranges::views::transform(
+            [](const FiniteElement& elem) {
+                return elem.nodes_id;
+            }),
+                std::ostream_iterator<std::vector<int>>(file));
+}
-- 
GitLab