Skip to content
Snippets Groups Projects
Commit b84b759a authored by Abdelkader Bouguettaiaa's avatar Abdelkader Bouguettaiaa
Browse files

avant duplication

parent c29157e4
No related branches found
No related tags found
No related merge requests found
Pipeline #2024 failed
Stage: test
Stage: coverage
Show whitespace changes
Inline Side-by-side
src/scheme/DDFVSolver.cpp
+ 100
− 14
View file @ b84b759a
......@@ -87,6 +87,7 @@ template <MeshConcept MeshType>
class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolver
{
private:
static constexpr size_t Dimension = MeshType::Dimension;
using Rdxd = TinyMatrix<Dimension>;
......@@ -174,15 +175,17 @@ class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolve
return bc_list;
}
/*
void _applySymmetryBCvelocity(const BoundaryConditionList& bc_list,
const MeshType& mesh,
const double& dt_moyen,
NodeValue<CellId>& computed_node_dual_cell_id,
CellValue<Rd>& u_r,
CellValue<double>& p_j,
CellValue<double>& q_j,
CellValue<double>& m_r
) const;
*/
void _applyDirichletBCvelocity(const BoundaryConditionList& bc_list,
NodeValue<CellId>& computed_node_dual_cell_id,
CellValue<Rd>& new_u_r,
......@@ -192,7 +195,7 @@ class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolve
const MeshType& mesh,
const double& dt_moyen,
NodeValue<CellId>& computed_node_dual_cell_id,
CellValue<Rd>& u_r,
CellValue<Rd>& flux_vitesse_dual,
CellValue<double>& p_j,
CellValue<double>& q_j,
CellValue<double>& m_r
......@@ -214,8 +217,15 @@ class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolve
{
this->_applyDirichletBCvelocity(bc_list, computed_node_dual_cell_id, new_u_r, flux_vitesse_dual);
this->_applyPeriodicBCvelocity(bc_list, mesh, dt_moyen, computed_node_dual_cell_id, flux_vitesse_dual, p_j, q_j, m_r);
//this->_applySymmetryBCvelocity(bc_list, mesh, dt_moyen, computed_node_dual_cell_id, flux_vitesse_dual, p_j, q_j, m_r);
}
/*
void _applySymmetryBCenergie(const BoundaryConditionList& bc_list,
const MeshType& mesh_p,
NodeValue<CellId>& computed_node_dual_cell_id,
CellValue<double>& flux_energie_dual,
CellValue<Rd>& u_r) const;
*/
void _applyDirichletBCenergie(const BoundaryConditionList& bc_list,
const MeshType& mesh_p,
NodeValue<CellId>& computed_node_dual_cell_id,
......@@ -238,6 +248,7 @@ class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolve
{
this->_applyDirichletBCenergie(bc_list, mesh, computed_node_dual_cell_id, flux_energie_dual, u_r);
this->_applyPeriodicBCenergie(bc_list, mesh, computed_node_dual_cell_id, flux_energie_dual, u_j);
}
void _applyPeriodicBCdensite(const BoundaryConditionList& bc_list,
......@@ -289,6 +300,7 @@ class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolve
const DiscreteScalarFunction& q_r,
const BoundaryConditionList& bc_list) const
{
//const auto& node_to_face_matrix_p = mesh_p.connectivity().nodeToFaceMatrix();
TinyMatrix<Dimension> I =identity;
double gamma = 1.4;
......@@ -640,9 +652,9 @@ class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolve
//double lambda=new_mesh_data_d.sumOverRLjr()[dual_cell_id_r];
//double lambda = new_Vr[dual_cell_id_r];
if(div_u_r[dual_cell_id_r]>0.0){
new_q_r[dual_cell_id_r] = - 0.0 * new_Vr[dual_cell_id_r] * new_rho_r[dual_cell_id_r] * div_u_r[dual_cell_id_r] * sqrt(gamma * new_p_r[dual_cell_id_r] / new_rho_r[dual_cell_id_r]);
new_q_r[dual_cell_id_r] = - 1.0 * new_Vr[dual_cell_id_r] * new_rho_r[dual_cell_id_r] * div_u_r[dual_cell_id_r] * sqrt(gamma * new_p_r[dual_cell_id_r] / new_rho_r[dual_cell_id_r]);
}else{
new_q_r[dual_cell_id_r] = - 0.0 * new_Vr[dual_cell_id_r] * new_rho_r[dual_cell_id_r] * div_u_r[dual_cell_id_r] * sqrt(gamma * new_p_r[dual_cell_id_r] / new_rho_r[dual_cell_id_r]);
new_q_r[dual_cell_id_r] = - 1.0 * new_Vr[dual_cell_id_r] * new_rho_r[dual_cell_id_r] * div_u_r[dual_cell_id_r] * sqrt(gamma * new_p_r[dual_cell_id_r] / new_rho_r[dual_cell_id_r]);
}
//new_q_r[dual_cell_id_r] = -lambda * new_rho_r[dual_cell_id_r] * (-0.0 * c + lambda * grad_u_r[dual_cell_id_r]) * grad_u_r[dual_cell_id_r];
};
......@@ -666,9 +678,9 @@ class DDFVSolverHandler::DDFVSolver final : public DDFVSolverHandler::IDDFVSolve
//double lambda = new_mesh_data_p.sumOverRLjr()[j];
//double lambda = new_Vj[j];
if (div_u_j[j]>0.0){
new_q_j[j] = - 0.0 * new_Vj[j] * new_rho_j[j] * div_u_j[j] * sqrt(gamma * new_p_j[j] / new_rho_j[j]);
new_q_j[j] = - 1.0 * new_Vj[j] * new_rho_j[j] * div_u_j[j] * sqrt(gamma * new_p_j[j] / new_rho_j[j]);
}else{
new_q_j[j] = -0.0 * new_Vj[j] * new_rho_j[j] * div_u_j[j] * sqrt(gamma * new_p_j[j] / new_rho_j[j]);
new_q_j[j] = -1.0 * new_Vj[j] * new_rho_j[j] * div_u_j[j] * sqrt(gamma * new_p_j[j] / new_rho_j[j]);
}
//new_q_j[j] = -lambda * new_rho_j[j] * (-0.0 * c + lambda * grad_u_j[j]) * grad_u_j[j];
......@@ -1032,6 +1044,75 @@ DDFVSolverHandler::DDFVSolver<MeshType>::_applyDirichletBCvelocity(const Boundar
}
}
/*
template <MeshConcept MeshType>
void
DDFVSolverHandler::DDFVSolver<MeshType>::_applySymmetryBCenergie(const BoundaryConditionList& bc_list,
const MeshType& mesh_p,
NodeValue<CellId>& computed_node_dual_cell_id,
CellValue<double>& flux_energie_dual,
CellValue<Rd>& u_j) const
{
MeshDataType& mesh_data_p = MeshDataManager::instance().getMeshData(mesh_p);
const FaceValue<const TinyVector<Dimension>> primal_nl = mesh_data_p.nl();
const NodeValuePerCell<const Rd> Cjr = mesh_data_p.Cjr();
const auto& node_to_cell_matrix_p = mesh_p.connectivity().nodeToCellMatrix();
const auto node_local_numbers_in_their_cells = mesh_p.connectivity().nodeLocalNumbersInTheirCells();
Rd flux = zero;
double masse = 0;
if(Dimension==2){
const auto& node_to_face_matrix_p = mesh_p.connectivity().nodeToFaceMatrix();
for (const auto& boundary_condition : bc_list) {
std::visit(
[&](auto&& bc) {
using T = std::decay_t<decltype(bc)>;
if constexpr (std::is_same_v<SymmetryBoundaryCondition, T>) {
const auto& node_list = bc.nodeList();
//std::cout<<"incroyable ça marche oidgidfxdizdjqpizjdzizqdijoqzlkd"<<std::endl;
for(size_t i_node = 0 ; i_node < node_list.size() ; i_node++ ){
NodeId node_id = node_list[i_node];
const auto& node_faces = node_to_face_matrix_p[node_id];
const auto& node_cells = node_to_cell_matrix_p[node_id];
CellId node_dual_cell_id = computed_node_dual_cell_id[node_id];
const auto& node_local_number_in_its_cells = node_local_numbers_in_their_cells.itemArray(node_id);
//flux_vitesse_dual[node_dual_cell_id] = zero;
for(size_t f = 0; f < node_faces.size(); f++){
FaceId face_id = node_faces[f];
if(mesh_p.connectivity().isBoundaryFace()[face_id]==true){
Rd face_b;
face_b[0]=-primal_nl[face_id][1];
face_b[1]=primal_nl[face_id][0];
for(size_t j = 0; j < node_cells.size(); j++){
CellId cell_id = node_cells[j];
const unsigned int R = node_local_number_in_its_cells[j];
//flux_energie_dual[node_dual_cell_id] += p_j[cell_id] * dot(Cjr(j,R), face_b) * face_b;
}
}
//std::cout<<node_id<<"et à coté la maille associé"<<j<<"et la numérotation local"<< R <<"et pour la forme le Cjr "<<flux<<std::endl;
//std::cout<<node_cells.size()<<"=1"<<std::endl;
}
}
}
},
boundary_condition);
}
}
}
template <MeshConcept MeshType>
void
DDFVSolverHandler::DDFVSolver<MeshType>::_applySymmetryBCvelocity(const BoundaryConditionList& bc_list,
......@@ -1040,6 +1121,7 @@ DDFVSolverHandler::DDFVSolver<MeshType>::_applySymmetryBCvelocity(const Boundary
NodeValue<CellId>& computed_node_dual_cell_id,
CellValue<Rd>& flux_vitesse_dual,
CellValue<double>& p_j,
CellValue<double>& q_j,
CellValue<double>& m_r) const
{
MeshDataType& mesh_data_p = MeshDataManager::instance().getMeshData(mesh_p);
......@@ -1049,8 +1131,10 @@ DDFVSolverHandler::DDFVSolver<MeshType>::_applySymmetryBCvelocity(const Boundary
const auto node_local_numbers_in_their_cells = mesh_p.connectivity().nodeLocalNumbersInTheirCells();
Rd flux = zero;
double masse = 0;
if(Dimension==2){
const auto& node_to_face_matrix_p = mesh_p.connectivity().nodeToFaceMatrix();
//l'idée est de bouclé sur les segments ou faces autours du sommets et de vérifier si il est sur le bords si il l'est c'est bon on fait tout suite à le faire 2 fois il y a une multiplication par 2 de toute facon qui doit apparaitre
......@@ -1073,17 +1157,17 @@ DDFVSolverHandler::DDFVSolver<MeshType>::_applySymmetryBCvelocity(const Boundary
flux_vitesse_dual[node_dual_cell_id] = zero;
for(size_t f = 0; f < node_faces.size(); f++)
if(mesh_p.connectivity().isBoundaryFace()[f]==true){
for(size_t f = 0; f < node_faces.size(); f++){
FaceId face_id = node_faces[f];
if(mesh_p.connectivity().isBoundaryFace()[face_id]==true){
Rd face_b;
face_b[0]=-primal_nl[face_id][1];
face_b[1]=primal_nl[face_id][0];
for(size_t j = 0; j < node_cells.size(); j++){
CellId cell_id = node_cells[j];
const unsigned int R = node_local_number_in_its_cells[j];
flux_vitesse_dual[node_dual_cell_id] += dt_moyen * p_j[cell_id] * dot(Cjr(j,R), face_b) * face_b;
flux_vitesse_dual[node_dual_cell_id] += dt_moyen * (p_j[cell_id] + q_j[cell_id]) * dot(Cjr(j,R), face_b) * face_b;
}
}
......@@ -1092,11 +1176,13 @@ DDFVSolverHandler::DDFVSolver<MeshType>::_applySymmetryBCvelocity(const Boundary
//std::cout<<node_cells.size()<<"=1"<<std::endl;
}
}
}
},
boundary_condition);
}
}
}
*/
template <MeshConcept MeshType>
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment