From d432e9d09e4e4db160ad96930e2ee60d1e97c265 Mon Sep 17 00:00:00 2001 From: Drouard <axelle.drouard@orange.fr> Date: Thu, 21 Jul 2022 14:05:58 +0200 Subject: [PATCH] Multiply the scheme by Vj --- src/scheme/ScalarNodalScheme.cpp | 42 +++++++++++++++----------------- 1 file changed, 19 insertions(+), 23 deletions(-) diff --git a/src/scheme/ScalarNodalScheme.cpp b/src/scheme/ScalarNodalScheme.cpp index 25bce53f9..f9767abd0 100644 --- a/src/scheme/ScalarNodalScheme.cpp +++ b/src/scheme/ScalarNodalScheme.cpp @@ -761,10 +761,10 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand const size_t i_node_k = node_local_number_in_its_cells[i_cell_k]; S1(cell_id_j, cell_id_k) += - dt * (1. - cn_coeff / 2.) / Vj[cell_id_j] * + dt * (1. - cn_coeff / 2.) * dot(kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); S2(cell_id_j, cell_id_k) -= - dt * (cn_coeff / 2.) / Vj[cell_id_j] * + dt * (cn_coeff / 2.) * dot(kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); } } @@ -783,12 +783,12 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand const size_t i_node_k = node_local_number_in_its_cells[i_cell_k]; S1(cell_id_j, cell_id_k) += - dt * (1. - cn_coeff / 2.) / Vj[cell_id_j] * + dt * (1. - cn_coeff / 2.) * dot(kappar_invBetar[node_id] * (Cjr(cell_id_k, i_node_k) - theta(cell_id_k, i_node_k) / sum_theta[node_id] * sum_Cjr[node_id]), Q * Cjr(cell_id_j, i_node_j)); S2(cell_id_j, cell_id_k) -= - dt * (cn_coeff / 2.) / Vj[cell_id_j] * + dt * (cn_coeff / 2.) * dot(kappar_invBetar[node_id] * (Cjr(cell_id_k, i_node_k) - theta(cell_id_k, i_node_k) / sum_theta[node_id] * sum_Cjr[node_id]), Q * Cjr(cell_id_j, i_node_j)); @@ -804,10 +804,10 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand const size_t i_node_k = node_local_number_in_its_cells[i_cell_k]; S1(cell_id_j, cell_id_k) += - dt * (1. - cn_coeff / 2.) / Vj[cell_id_j] * + dt * (1. - cn_coeff / 2.) * dot(kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); S2(cell_id_j, cell_id_k) -= - dt * (cn_coeff / 2.) / Vj[cell_id_j] * + dt * (cn_coeff / 2.) * dot(kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); } } @@ -821,10 +821,10 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand const size_t i_node_k = node_local_number_in_its_cells[i_cell_k]; S1(cell_id_j, cell_id_k) += - dt * (1. - cn_coeff / 2.) / Vj[cell_id_j] * + dt * (1. - cn_coeff / 2.) * dot(corner_kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); S2(cell_id_j, cell_id_k) -= - dt * (cn_coeff / 2.) / Vj[cell_id_j] * + dt * (cn_coeff / 2.) * dot(corner_kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); } } @@ -832,8 +832,8 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand } for (CellId cell_id = 0; cell_id < mesh->numberOfCells(); ++cell_id) { - S1(cell_id, cell_id) += 1; - S2(cell_id, cell_id) += 1; + S1(cell_id, cell_id) += Vj[cell_id]; + S2(cell_id, cell_id) += Vj[cell_id]; } CellValue<const double> fj = f->cellValues(); @@ -852,13 +852,13 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand Vector<double> b{mesh->numberOfCells()}; b = zero; for (CellId cell_id = 0; cell_id < mesh->numberOfCells(); ++cell_id) { - b[cell_id] = fj[cell_id]; + b[cell_id] = Vj[cell_id] * fj[cell_id]; }; Vector<double> b_prev{mesh->numberOfCells()}; b_prev = zero; for (CellId cell_id = 0; cell_id < mesh->numberOfCells(); ++cell_id) { - b_prev[cell_id] = fj_prev[cell_id]; + b_prev[cell_id] = Vj[cell_id] * fj_prev[cell_id]; }; for (NodeId node_id = 0; node_id < mesh->numberOfNodes(); ++node_id) { @@ -874,11 +874,11 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand const CellId cell_id = node_to_cell[i_cell]; const size_t i_node = node_local_number_in_its_cells[i_cell]; - b[cell_id] += 1. / Vj[cell_id] * node_boundary_values[node_id] * + b[cell_id] += node_boundary_values[node_id] * (1. / (sum_theta[node_id] * l2Norm(node_kappar[node_id] * n)) * dot(kappar_invBetar[node_id] * sum_Cjr[node_id], Q * Cjr(cell_id, i_node)) + dot(Cjr(cell_id, i_node), n)); - b_prev[cell_id] += 1. / Vj[cell_id] * node_boundary_prev_values[node_id] * + b_prev[cell_id] += node_boundary_prev_values[node_id] * (1. / (sum_theta[node_id] * l2Norm(node_kappar[node_id] * n)) * dot(kappar_invBetar[node_id] * sum_Cjr[node_id], Q * Cjr(cell_id, i_node)) + dot(Cjr(cell_id, i_node), n)); @@ -891,8 +891,8 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand FaceId face_id = node_to_face[i_face]; sum_mes_l += mes_l[face_id]; } - b[cell_id] += 0.5 / Vj[cell_id] * node_boundary_values[node_id] * sum_mes_l; - b_prev[cell_id] += 0.5 / Vj[cell_id] * node_boundary_prev_values[node_id] * sum_mes_l; + b[cell_id] += 0.5 * node_boundary_values[node_id] * sum_mes_l; + b_prev[cell_id] += 0.5 * node_boundary_prev_values[node_id] * sum_mes_l; } } else if (node_is_dirichlet[node_id]) { if (not node_is_corner[node_id]) { @@ -904,11 +904,9 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand const CellId cell_id_k = node_to_cell[i_cell_k]; const size_t i_node_k = node_local_number_in_its_cells[i_cell_k]; - b[cell_id_j] += 1. / Vj[cell_id_j] * - dot(node_boundary_values[node_id] * kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), + b[cell_id_j] += dot(node_boundary_values[node_id] * kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); b_prev[cell_id_j] += - 1. / Vj[cell_id_j] * dot(node_boundary_prev_values[node_id] * kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); } @@ -924,11 +922,9 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand const size_t i_node_k = node_local_number_in_its_cells[i_cell_k]; b[cell_id_j] += - 1. / Vj[cell_id_j] * dot(node_boundary_values[node_id] * corner_kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); b_prev[cell_id_j] += - 1. / Vj[cell_id_j] * dot(node_boundary_prev_values[node_id] * corner_kappar_invBetar[node_id] * Cjr(cell_id_k, i_node_k), Cjr(cell_id_j, i_node_j)); } @@ -963,7 +959,7 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand ur[node_id] += Pj[cell_id] * Cjr(cell_id, i_node); } ur[node_id] = -inverse(node_betar[node_id]) * ur[node_id]; - std::cout << node_id << "; ur = " << ur[node_id] << "\n"; + // std::cout << node_id << "; ur = " << ur[node_id] << "\n"; } else if (is_boundary_node[node_id] and node_is_dirichlet[node_id]) { for (size_t i_cell = 0; i_cell < node_to_cell.size(); ++i_cell) { const CellId cell_id = node_to_cell[i_cell]; @@ -975,7 +971,7 @@ class ScalarNodalSchemeHandler::ScalarNodalScheme : public ScalarNodalSchemeHand } else { ur[node_id] = inverse(corner_betar[node_id]) * ur[node_id]; } - std::cout << "bord : " << node_id << "; ur = " << ur[node_id] << "\n"; + // std::cout << "bord : " << node_id << "; ur = " << ur[node_id] << "\n"; } } -- GitLab