From 1151fb33b26dcc545aa0a4e26961cea9ddaa2252 Mon Sep 17 00:00:00 2001 From: "dino.celebic" Date: Mon, 26 Jan 2026 17:56:16 +0100 Subject: [PATCH] made better --- mgrid_2/geom.cpp | 38 +++++++++++++++++++++++++++++++------- mgrid_2/geom.h | 4 +++- mgrid_2/main.cpp | 11 +++++++---- solid-cpp/ownSolver.cpp | 17 +++++++++++++++-- 4 files changed, 56 insertions(+), 14 deletions(-) diff --git a/mgrid_2/geom.cpp b/mgrid_2/geom.cpp index 254efc3..08bf123 100644 --- a/mgrid_2/geom.cpp +++ b/mgrid_2/geom.cpp @@ -496,14 +496,13 @@ void Mesh::Visualize_paraview(vector const &v) const return; } -std::tuple Mesh::AverageVectorFunction_perSubdomain(std::vector &v, int target_sd) const +double Mesh::AverageVectorFunction_perSubdomain(std::vector &v, int target_sd) const { assert(2==Ndims()); int const nnode = Nnodes(); // number of vertices in mesh assert( nnode == static_cast(v.size()) ); double cumulative_element_temp = 0.0; - double elements_temp_reached = 0.0; int subdomain_element_counter = 0; for (int e = 0; e < Nelems(); ++e) // loop over all elements { @@ -519,16 +518,41 @@ std::tuple Mesh::AverageVectorFunction_perSubdomain(std::vector 67.0) // check if element has reached temp + return cumulative_element_temp/subdomain_element_counter; +} + +double Mesh::CheckTemp_mult(std::vector &v, int target_sd, double goal_temp) const +{ + assert(2==Ndims()); + int const nnode = Nnodes(); // number of vertices in mesh + assert( nnode == static_cast(v.size()) ); + + double elements_temp_reached = 0.0; + int subdomain_element_counter = 0; + for (int e = 0; e < Nelems(); ++e) // loop over all elements + { + int sd = ElementSubdomains[e]; // get subdomain of element e + if (sd == target_sd) // if is target subdomain then + { + subdomain_element_counter++; + int base = e * _nvert_e; // get starting index of element in coordinate vector + double cumulative_node_temp = 0.0; + for (int k = 0; k < _nvert_e; ++k) // loop over vertices of element { - elements_temp_reached++; // add to counter + int node = _ia[base + k]; // global index of vertex + cumulative_node_temp += v[node]; // set function + } + if (cumulative_node_temp/_nvert_e > goal_temp) + { + elements_temp_reached++; } } } - // average temperature % of elements reached temperature 67 - return std::make_tuple(cumulative_element_temp/subdomain_element_counter, 100 * elements_temp_reached / static_cast(subdomain_element_counter)); + + return 100 * elements_temp_reached / subdomain_element_counter; } diff --git a/mgrid_2/geom.h b/mgrid_2/geom.h index d32880e..84717c5 100644 --- a/mgrid_2/geom.h +++ b/mgrid_2/geom.h @@ -231,8 +231,10 @@ public: [[nodiscard]] virtual std::vector Index_DirichletNodes() const; - [[nodiscard]] std::tuple AverageVectorFunction_perSubdomain(std::vector &v, int target_sd) const; + [[nodiscard]] double AverageVectorFunction_perSubdomain(std::vector &v, int target_sd) const; + [[nodiscard]] double CheckTemp_mult(std::vector &v, int target_sd, double goal_temp) const; + /** * Determines the indices of those vertices with Dirichlet boundary conditions. * diff --git a/mgrid_2/main.cpp b/mgrid_2/main.cpp index 8f4037e..ea5f399 100644 --- a/mgrid_2/main.cpp +++ b/mgrid_2/main.cpp @@ -84,10 +84,12 @@ int main(int argc, char **argv ) auto t3 = system_clock::now(); // start timer double average_cup_temperature = u0_mug; double percentage_temp_reached = 0.0; + double goal_temp = 67.0; + double goal_perc = 60.0; double time_count = 0; - while (average_cup_temperature < 67.0) - // while (percentage_temp_reached < 60.0) + while (average_cup_temperature < goal_temp) + // while (percentage_temp_reached < goal_perc) //for (int step = 0; step < steps; ++step) { vector G(Mdt.Nrows(), 0.0); @@ -102,9 +104,10 @@ int main(int argc, char **argv ) JacobiSolve(SK, H, uv); // solve: (M/dt + K + C) * u_{n+1} = F + M/dt * u_{n} // ----- SK ----- ------ H ------- - tie(average_cup_temperature, percentage_temp_reached) = mesh_c.AverageVectorFunction_perSubdomain(uv, 0); + average_cup_temperature = mesh_c.AverageVectorFunction_perSubdomain(uv, 0); + percentage_temp_reached = mesh_c.CheckTemp_mult(uv, 0, goal_temp); cout << "Average cup temperature: " << average_cup_temperature << " after " << time_count << " seconds. " << endl; - cout << "% of mug elements reached temperature 67º: " << percentage_temp_reached << endl; + cout << "% of elements reached temperature " << goal_temp << "ºC: " << percentage_temp_reached << endl; time_count += dt; } auto t4 = system_clock::now(); // stop timer diff --git a/solid-cpp/ownSolver.cpp b/solid-cpp/ownSolver.cpp index 6ecfe78..c0c93e3 100644 --- a/solid-cpp/ownSolver.cpp +++ b/solid-cpp/ownSolver.cpp @@ -88,9 +88,13 @@ int main(int argc, char **argv ) auto t3 = system_clock::now(); // start timer double average_cup_temperature = u0_mug; + double percentage_temp_reached = 0.0; + double goal_temp = 67.0; + double goal_perc = 60.0; double time_count = 0; - while (average_cup_temperature < 67.0) + while (average_cup_temperature < goal_temp) + // while (percentage_temp_reached < goal_perc) //for (int step = 0; step < steps; ++step) { vector G(Mdt.Nrows(), 0.0); @@ -106,7 +110,9 @@ int main(int argc, char **argv ) // ----- SK ----- ------ H ------- average_cup_temperature = mesh_c.AverageVectorFunction_perSubdomain(uv, 0); + percentage_temp_reached = mesh_c.CheckTemp_mult(uv, 0, goal_temp); cout << "Average cup temperature: " << average_cup_temperature << " after " << time_count << " seconds. " << endl; + cout << "% of elements reached temperature " << goal_temp << "ºC: " << percentage_temp_reached << endl; time_count += dt; } auto t4 = system_clock::now(); // stop timer @@ -124,6 +130,7 @@ int main(int argc, char **argv ) t3 = system_clock::now(); // start timer double average_coffee_temperature = u0_coffee; + percentage_temp_reached = 0.0; // ------------------------ initialize preCICE ------------------------ int commRank = 0; @@ -162,8 +169,12 @@ int main(int argc, char **argv ) // ------------------------ timestepping ------------------------ + goal_temp = 50.0; + goal_perc = 60.0; + time_count = 0; - while (average_coffee_temperature > 50.0) + while (average_cup_temperature > goal_temp) + // while (percentage_temp_reached > goal_perc) { preciceDt = participantSolid.getMaxTimeStepSize(); solverDt = 1.0; @@ -190,7 +201,9 @@ int main(int argc, char **argv ) // ----- SK ----- ------ H ------- average_coffee_temperature = mesh_c.AverageVectorFunction_perSubdomain(uv, 1); + percentage_temp_reached = mesh_c.CheckTemp_mult(uv, 0, goal_temp); cout << "Average coffee temperature: " << average_coffee_temperature << " after " << time_count << " seconds. " << endl; + cout << "% of elements reached temperature " << goal_temp << "ºC: " << percentage_temp_reached << endl; // ----- write the heat-flux, so openFOAM can read it {