From 1d4ee11edbdd1e79134160fc4428f6d8f03901a4 Mon Sep 17 00:00:00 2001
From: HOCH PHILIPPE <philippe.hoch@gmail.com>
Date: Tue, 1 Apr 2025 10:08:57 +0200
Subject: [PATCH] New other implementation for Wall boundary at edge in 3D
---
...eViscousFormEulerianCompositeSolver_v2.cpp | 47 ++++++++++----
...scousFormEulerianCompositeSolver_v2_o2.cpp | 46 ++++++++++----
src/scheme/RusanovEulerianCompositeSolver.cpp | 62 ++++++++++++++-----
.../RusanovEulerianCompositeSolver_o2.cpp | 47 ++++++++++----
.../RusanovEulerianCompositeSolver_v2.cpp | 48 ++++++++++----
.../RusanovEulerianCompositeSolver_v2_o2.cpp | 48 ++++++++++----
6 files changed, 219 insertions(+), 79 deletions(-)
diff --git a/src/scheme/RoeViscousFormEulerianCompositeSolver_v2.cpp b/src/scheme/RoeViscousFormEulerianCompositeSolver_v2.cpp
index 397eb85e7..fd856ffbc 100644
--- a/src/scheme/RoeViscousFormEulerianCompositeSolver_v2.cpp
+++ b/src/scheme/RoeViscousFormEulerianCompositeSolver_v2.cpp
@@ -620,26 +620,49 @@ class RoeViscousFormEulerianCompositeSolver_v2
const auto& edge_local_numbers_in_their_cells = mesh.connectivity().edgeLocalNumbersInTheirCells();
- const auto& edge_list = bc.edgeList();
+ const auto& edge_to_face_matrix = mesh.connectivity().edgeToFaceMatrix();
+
+ const auto& edge_local_numbers_in_their_faces = mesh.connectivity().edgeLocalNumbersInTheirFaces();
- const auto Cje = mesh_data.Cje();
+ const auto& edge_list = bc.edgeList();
for (size_t i_edge = 0; i_edge < edge_list.size(); ++i_edge) {
const EdgeId edge_id = edge_list[i_edge];
- const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
-
- const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
+ const auto& edge_face_list = edge_to_face_matrix[edge_id];
- // Normal locale approchée
+ // on va chercher les normale d'une face issue du edge de CL et contenue dans le faceList
Rd normal(zero);
- for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
- CellId edge_cell_id = edge_cell_list[i_cell];
- size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
- normal += Cje(edge_cell_id, edge_local_number_in_cell);
+ int nbnormal = 0;
+ for (size_t i_face = 0; i_face < edge_face_list.size(); ++i_face) {
+ FaceId edge_face_id = edge_face_list[i_face];
+
+ for (size_t i_facebc = 0; i_facebc < face_list.size(); ++i_facebc) {
+ const FaceId facebc_id = face_list[i_facebc];
+ if (edge_face_id == facebc_id) {
+ const auto& face_cell_list = face_to_cell_matrix[facebc_id];
+ Assert(face_cell_list.size() == 1);
+
+ CellId face_cell_id = face_cell_list[0];
+ size_t face_local_number_in_cell = face_local_numbers_in_their_cells(facebc_id, 0);
+
+ // Normal locale approchée
+ Rd normalloc = Cjf(face_cell_id, face_local_number_in_cell);
+ normalloc *= 1. / l2Norm(normalloc);
+ normal += normalloc;
+ ++nbnormal;
+ break;
+ }
+ }
}
- normal *= 1. / edge_cell_list.size();
- normal *= 1. / l2Norm(normal);
+
+ if (nbnormal == 0)
+ continue;
+ normal *= 1. / nbnormal;
+
+ const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
+
+ const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
CellId edge_cell_id = edge_cell_list[i_cell];
diff --git a/src/scheme/RoeViscousFormEulerianCompositeSolver_v2_o2.cpp b/src/scheme/RoeViscousFormEulerianCompositeSolver_v2_o2.cpp
index 2c9e1711c..934c39b96 100644
--- a/src/scheme/RoeViscousFormEulerianCompositeSolver_v2_o2.cpp
+++ b/src/scheme/RoeViscousFormEulerianCompositeSolver_v2_o2.cpp
@@ -625,26 +625,48 @@ class RoeViscousFormEulerianCompositeSolver_v2_o2
const auto& edge_local_numbers_in_their_cells = mesh.connectivity().edgeLocalNumbersInTheirCells();
- const auto& edge_list = bc.edgeList();
+ const auto& edge_to_face_matrix = mesh.connectivity().edgeToFaceMatrix();
+
+ const auto& edge_local_numbers_in_their_faces = mesh.connectivity().edgeLocalNumbersInTheirFaces();
- const auto Cje = mesh_data.Cje();
+ const auto& edge_list = bc.edgeList();
for (size_t i_edge = 0; i_edge < edge_list.size(); ++i_edge) {
const EdgeId edge_id = edge_list[i_edge];
- const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
+ const auto& edge_face_list = edge_to_face_matrix[edge_id];
- const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
-
- // Normal locale approchée
+ // on va chercher les normale d'une face issue du edge de CL et contenue dans le faceList
Rd normal(zero);
- for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
- CellId edge_cell_id = edge_cell_list[i_cell];
- size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
- normal += Cje(edge_cell_id, edge_local_number_in_cell);
+ int nbnormal = 0;
+ for (size_t i_face = 0; i_face < edge_face_list.size(); ++i_face) {
+ FaceId edge_face_id = edge_face_list[i_face];
+
+ for (size_t i_facebc = 0; i_facebc < face_list.size(); ++i_facebc) {
+ const FaceId facebc_id = face_list[i_facebc];
+ if (edge_face_id == facebc_id) {
+ const auto& face_cell_list = face_to_cell_matrix[facebc_id];
+ Assert(face_cell_list.size() == 1);
+
+ CellId face_cell_id = face_cell_list[0];
+ size_t face_local_number_in_cell = face_local_numbers_in_their_cells(facebc_id, 0);
+
+ // Normal locale approchée
+ Rd normalloc = Cjf(face_cell_id, face_local_number_in_cell);
+ normalloc *= 1. / l2Norm(normalloc);
+ normal += normalloc;
+ ++nbnormal;
+ break;
+ }
+ }
}
- normal *= 1. / edge_cell_list.size();
- normal *= 1. / l2Norm(normal);
+
+ if (nbnormal == 0)
+ continue;
+ normal *= 1. / nbnormal;
+ const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
+
+ const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
CellId edge_cell_id = edge_cell_list[i_cell];
diff --git a/src/scheme/RusanovEulerianCompositeSolver.cpp b/src/scheme/RusanovEulerianCompositeSolver.cpp
index 2bedc8879..8374ae716 100644
--- a/src/scheme/RusanovEulerianCompositeSolver.cpp
+++ b/src/scheme/RusanovEulerianCompositeSolver.cpp
@@ -285,27 +285,51 @@ class RusanovEulerianCompositeSolver
const auto& edge_local_numbers_in_their_cells = mesh.connectivity().edgeLocalNumbersInTheirCells();
- const auto& edge_list = bc.edgeList();
-
- const auto Cje = mesh_data.Cje();
+ const auto& edge_to_face_matrix = mesh.connectivity().edgeToFaceMatrix();
- for (size_t i_edge = 0; i_edge < edge_list.size(); ++i_edge) {
- const EdgeId edge_id = edge_list[i_edge];
+ const auto& edge_local_numbers_in_their_faces = mesh.connectivity().edgeLocalNumbersInTheirFaces();
- const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
+ const auto& edge_list = bc.edgeList();
- const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
+ for (size_t i_edge = 0; i_edge < edge_list.size(); ++i_edge) {
+ const EdgeId edge_id = edge_list[i_edge];
+ const auto& edge_face_list = edge_to_face_matrix[edge_id];
- // Normal locale approchée
+ // on va chercher les normale d'une face issue du edge de CL et contenue dans le faceList
Rd normal(zero);
- for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
- CellId edge_cell_id = edge_cell_list[i_cell];
- size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
- normal += Cje(edge_cell_id, edge_local_number_in_cell);
+ int nbnormal = 0;
+ for (size_t i_face = 0; i_face < edge_face_list.size(); ++i_face) {
+ FaceId edge_face_id = edge_face_list[i_face];
+
+ for (size_t i_facebc = 0; i_facebc < face_list.size(); ++i_facebc) {
+ const FaceId facebc_id = face_list[i_facebc];
+ if (edge_face_id == facebc_id) {
+ const auto& face_cell_list = face_to_cell_matrix[facebc_id];
+ Assert(face_cell_list.size() == 1);
+
+ CellId face_cell_id = face_cell_list[0];
+ size_t face_local_number_in_cell = face_local_numbers_in_their_cells(facebc_id, 0);
+
+ // Normal locale approchée
+ Rd normalloc = Cjf(face_cell_id, face_local_number_in_cell);
+ normalloc *= 1. / l2Norm(normalloc);
+ normal += normalloc;
+ ++nbnormal;
+ break;
+ }
+ }
}
- normal *= 1. / edge_cell_list.size();
+
+ if (nbnormal == 0)
+ continue;
+ normal *= 1. / nbnormal;
+
normal *= 1. / l2Norm(normal);
+ const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
+
+ const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
+
for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
CellId edge_cell_id = edge_cell_list[i_cell];
size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
@@ -381,7 +405,8 @@ class RusanovEulerianCompositeSolver
for (size_t dim = 0; dim < Dimension; ++dim)
stateNode[node_cell_id][node_local_number_in_cell][dim + 1] = vectorSym[dim];
- // stateNode[node_cell_id][node_local_number_in_cell][dim] = 0; // node_array_list[i_node][dim];
+ // stateNode[node_cell_id][node_local_number_in_cell][dim] = 0; //
+ node_array_list[i_node][dim];
}
}
@@ -493,7 +518,8 @@ class RusanovEulerianCompositeSolver
for (size_t dim = 0; dim < Dimension; ++dim)
stateNode[node_cell_id][node_local_number_in_cell][dim + 1] = vectorSym[dim];
- // stateNode[node_cell_id][node_local_number_in_cell][dim] = 0; // node_array_list[i_node][dim];
+ // stateNode[node_cell_id][node_local_number_in_cell][dim] = 0; //
+ // node_array_list[i_node][dim];
}
}
@@ -600,7 +626,8 @@ class RusanovEulerianCompositeSolver
for (size_t dim = 0; dim < Dimension; ++dim)
stateNode[node_cell_id][node_local_number_in_cell][dim + 1] = vectorSym[dim];
- // stateNode[node_cell_id][node_local_number_in_cell][dim] = 0; // node_array_list[i_node][dim];
+ // stateNode[node_cell_id][node_local_number_in_cell][dim] = 0; //
+ // node_array_list[i_node][dim];
}
}
@@ -1242,7 +1269,8 @@ class RusanovEulerianCompositeSolver
}
/*
std::cout << " maxErr face " << xl[l] << " err " << l2Norm(SumGjf) << " size face cell"
- << face_to_cell.size() << " Cell " << face_to_cell[0] << " / " << face_to_cell[1] << " Cjf0 "
+ << face_to_cell.size() << " Cell " << face_to_cell[0] << " / " << face_to_cell[1] << " Cjf0
+ "
<< Cjf(face_to_cell[0], 0) << " / Cjf1 " << Cjf(face_to_cell[1], 1) << "\n";
if (l2Norm(SumGjf) > 1e-3)
diff --git a/src/scheme/RusanovEulerianCompositeSolver_o2.cpp b/src/scheme/RusanovEulerianCompositeSolver_o2.cpp
index 6fbaafb7a..9e59e875d 100644
--- a/src/scheme/RusanovEulerianCompositeSolver_o2.cpp
+++ b/src/scheme/RusanovEulerianCompositeSolver_o2.cpp
@@ -293,26 +293,49 @@ class RusanovEulerianCompositeSolver_o2
const auto& edge_local_numbers_in_their_cells = mesh.connectivity().edgeLocalNumbersInTheirCells();
- const auto& edge_list = bc.edgeList();
+ const auto& edge_to_face_matrix = mesh.connectivity().edgeToFaceMatrix();
+
+ const auto& edge_local_numbers_in_their_faces = mesh.connectivity().edgeLocalNumbersInTheirFaces();
- const auto Cje = mesh_data.Cje();
+ const auto& edge_list = bc.edgeList();
for (size_t i_edge = 0; i_edge < edge_list.size(); ++i_edge) {
const EdgeId edge_id = edge_list[i_edge];
- const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
-
- const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
+ const auto& edge_face_list = edge_to_face_matrix[edge_id];
- // Normal locale approchée
+ // on va chercher les normale d'une face issue du edge de CL et contenue dans le faceList
Rd normal(zero);
- for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
- CellId edge_cell_id = edge_cell_list[i_cell];
- size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
- normal += Cje(edge_cell_id, edge_local_number_in_cell);
+ int nbnormal = 0;
+ for (size_t i_face = 0; i_face < edge_face_list.size(); ++i_face) {
+ FaceId edge_face_id = edge_face_list[i_face];
+
+ for (size_t i_facebc = 0; i_facebc < face_list.size(); ++i_facebc) {
+ const FaceId facebc_id = face_list[i_facebc];
+ if (edge_face_id == facebc_id) {
+ const auto& face_cell_list = face_to_cell_matrix[facebc_id];
+ Assert(face_cell_list.size() == 1);
+
+ CellId face_cell_id = face_cell_list[0];
+ size_t face_local_number_in_cell = face_local_numbers_in_their_cells(facebc_id, 0);
+
+ // Normal locale approchée
+ Rd normalloc = Cjf(face_cell_id, face_local_number_in_cell);
+ normalloc *= 1. / l2Norm(normalloc);
+ normal += normalloc;
+ ++nbnormal;
+ break;
+ }
+ }
}
- normal *= 1. / edge_cell_list.size();
- normal *= 1. / l2Norm(normal);
+
+ if (nbnormal == 0)
+ continue;
+ normal *= 1. / nbnormal;
+
+ const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
+
+ const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
CellId edge_cell_id = edge_cell_list[i_cell];
diff --git a/src/scheme/RusanovEulerianCompositeSolver_v2.cpp b/src/scheme/RusanovEulerianCompositeSolver_v2.cpp
index af7c70e37..52269034a 100644
--- a/src/scheme/RusanovEulerianCompositeSolver_v2.cpp
+++ b/src/scheme/RusanovEulerianCompositeSolver_v2.cpp
@@ -284,27 +284,49 @@ class RusanovEulerianCompositeSolver_v2
const auto& edge_local_numbers_in_their_cells = mesh.connectivity().edgeLocalNumbersInTheirCells();
- const auto& edge_list = bc.edgeList();
+ const auto& edge_to_face_matrix = mesh.connectivity().edgeToFaceMatrix();
+
+ const auto& edge_local_numbers_in_their_faces = mesh.connectivity().edgeLocalNumbersInTheirFaces();
- const auto Cje = mesh_data.Cje();
+ const auto& edge_list = bc.edgeList();
for (size_t i_edge = 0; i_edge < edge_list.size(); ++i_edge) {
- const EdgeId edge_id = edge_list[i_edge];
+ const EdgeId edge_id = edge_list[i_edge];
+ const auto& edge_face_list = edge_to_face_matrix[edge_id];
+
+ // on va chercher les normale d'une face issue du edge de CL et contenue dans le faceList
+ Rd normal(zero);
+ int nbnormal = 0;
+ for (size_t i_face = 0; i_face < edge_face_list.size(); ++i_face) {
+ FaceId edge_face_id = edge_face_list[i_face];
+
+ for (size_t i_facebc = 0; i_facebc < face_list.size(); ++i_facebc) {
+ const FaceId facebc_id = face_list[i_facebc];
+ if (edge_face_id == facebc_id) {
+ const auto& face_cell_list = face_to_cell_matrix[facebc_id];
+ Assert(face_cell_list.size() == 1);
+
+ CellId face_cell_id = face_cell_list[0];
+ size_t face_local_number_in_cell = face_local_numbers_in_their_cells(facebc_id, 0);
+
+ // Normal locale approchée
+ Rd normalloc = Cjf(face_cell_id, face_local_number_in_cell);
+ normalloc *= 1. / l2Norm(normalloc);
+ normal += normalloc;
+ ++nbnormal;
+ break;
+ }
+ }
+ }
+
+ if (nbnormal == 0)
+ continue;
+ normal *= 1. / nbnormal;
const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
- // Normal locale approchée .. need better def
- Rd normal(zero);
- for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
- CellId edge_cell_id = edge_cell_list[i_cell];
- size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
- normal += Cje(edge_cell_id, edge_local_number_in_cell);
- }
- normal *= 1. / edge_cell_list.size();
- normal *= 1. / l2Norm(normal);
-
for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
CellId edge_cell_id = edge_cell_list[i_cell];
size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
diff --git a/src/scheme/RusanovEulerianCompositeSolver_v2_o2.cpp b/src/scheme/RusanovEulerianCompositeSolver_v2_o2.cpp
index be4c0495e..399c0c230 100644
--- a/src/scheme/RusanovEulerianCompositeSolver_v2_o2.cpp
+++ b/src/scheme/RusanovEulerianCompositeSolver_v2_o2.cpp
@@ -292,27 +292,49 @@ class RusanovEulerianCompositeSolver_v2_o2
const auto& edge_local_numbers_in_their_cells = mesh.connectivity().edgeLocalNumbersInTheirCells();
- const auto& edge_list = bc.edgeList();
+ const auto& edge_to_face_matrix = mesh.connectivity().edgeToFaceMatrix();
+
+ const auto& edge_local_numbers_in_their_faces = mesh.connectivity().edgeLocalNumbersInTheirFaces();
- const auto Cje = mesh_data.Cje();
+ const auto& edge_list = bc.edgeList();
for (size_t i_edge = 0; i_edge < edge_list.size(); ++i_edge) {
- const EdgeId edge_id = edge_list[i_edge];
+ const EdgeId edge_id = edge_list[i_edge];
+ const auto& edge_face_list = edge_to_face_matrix[edge_id];
+
+ // on va chercher les normale d'une face issue du edge de CL et contenue dans le faceList
+ Rd normal(zero);
+ int nbnormal = 0;
+ for (size_t i_face = 0; i_face < edge_face_list.size(); ++i_face) {
+ FaceId edge_face_id = edge_face_list[i_face];
+
+ for (size_t i_facebc = 0; i_facebc < face_list.size(); ++i_facebc) {
+ const FaceId facebc_id = face_list[i_facebc];
+ if (edge_face_id == facebc_id) {
+ const auto& face_cell_list = face_to_cell_matrix[facebc_id];
+ Assert(face_cell_list.size() == 1);
+
+ CellId face_cell_id = face_cell_list[0];
+ size_t face_local_number_in_cell = face_local_numbers_in_their_cells(facebc_id, 0);
+
+ // Normal locale approchée
+ Rd normalloc = Cjf(face_cell_id, face_local_number_in_cell);
+ normalloc *= 1. / l2Norm(normalloc);
+ normal += normalloc;
+ ++nbnormal;
+ break;
+ }
+ }
+ }
+
+ if (nbnormal == 0)
+ continue;
+ normal *= 1. / nbnormal;
const auto& edge_cell_list = edge_to_cell_matrix[edge_id];
const auto& edge_local_number_in_its_cells = edge_local_numbers_in_their_cells.itemArray(edge_id);
- // Normal locale approchée
- Rd normal(zero);
- for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
- CellId edge_cell_id = edge_cell_list[i_cell];
- size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
- normal += Cje(edge_cell_id, edge_local_number_in_cell);
- }
- normal *= 1. / edge_cell_list.size();
- normal *= 1. / l2Norm(normal);
-
for (size_t i_cell = 0; i_cell < edge_cell_list.size(); ++i_cell) {
CellId edge_cell_id = edge_cell_list[i_cell];
size_t edge_local_number_in_cell = edge_local_number_in_its_cells[i_cell];
--
GitLab