From c39a976ebd995dde5b513a0ea66d76ed6f5a972c Mon Sep 17 00:00:00 2001 From: "dino.celebic" Date: Mon, 26 Jan 2026 17:24:33 +0100 Subject: [PATCH] add another 67 degree break condition --- generate_mesh/visualize_results.m | 1 + mgrid_2/geom.cpp | 14 ++++++++++---- mgrid_2/geom.h | 3 +-- mgrid_2/main.cpp | 5 ++++- mgrid_2/visualize_results.m | 1 + 5 files changed, 17 insertions(+), 7 deletions(-) diff --git a/generate_mesh/visualize_results.m b/generate_mesh/visualize_results.m index 0f3546b..9295c7a 100644 --- a/generate_mesh/visualize_results.m +++ b/generate_mesh/visualize_results.m @@ -23,6 +23,7 @@ h = patch('Faces', ia, 'Vertices', xc, 'FaceVertexCData', v, ... axis equal tight; colorbar; colormap(jet); +caxis([18 85]); title('Heat distribution'); waitfor(h) % wait for closing the figure \ No newline at end of file diff --git a/mgrid_2/geom.cpp b/mgrid_2/geom.cpp index 441778c..254efc3 100644 --- a/mgrid_2/geom.cpp +++ b/mgrid_2/geom.cpp @@ -496,13 +496,14 @@ void Mesh::Visualize_paraview(vector const &v) const return; } -double Mesh::AverageVectorFunction_perSubdomain(std::vector &v, int target_sd) const +std::tuple 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 { @@ -518,14 +519,19 @@ double Mesh::AverageVectorFunction_perSubdomain(std::vector &v, int targ cumulative_node_temp += v[node]; // set function } cumulative_element_temp += cumulative_node_temp/_nvert_e; + + + if (cumulative_node_temp/_nvert_e > 67.0) // check if element has reached temp + { + elements_temp_reached++; // add to counter + } } } - - return cumulative_element_temp/subdomain_element_counter; + // 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)); } - vector Mesh::Index_DirichletNodes() const { assert(2==Ndims()); // not in 3D currently diff --git a/mgrid_2/geom.h b/mgrid_2/geom.h index a7f0926..d32880e 100644 --- a/mgrid_2/geom.h +++ b/mgrid_2/geom.h @@ -231,9 +231,8 @@ public: [[nodiscard]] virtual std::vector Index_DirichletNodes() const; - [[nodiscard]] double AverageVectorFunction_perSubdomain(std::vector &v, int target_sd) const; + [[nodiscard]] std::tuple AverageVectorFunction_perSubdomain(std::vector &v, int target_sd) const; - /** * Determines the indices of those vertices with Dirichlet boundary conditions. * diff --git a/mgrid_2/main.cpp b/mgrid_2/main.cpp index 510a453..8f4037e 100644 --- a/mgrid_2/main.cpp +++ b/mgrid_2/main.cpp @@ -83,9 +83,11 @@ 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 time_count = 0; while (average_cup_temperature < 67.0) + // while (percentage_temp_reached < 60.0) //for (int step = 0; step < steps; ++step) { vector G(Mdt.Nrows(), 0.0); @@ -100,8 +102,9 @@ 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 ------- - average_cup_temperature = mesh_c.AverageVectorFunction_perSubdomain(uv, 0); + tie(average_cup_temperature, percentage_temp_reached) = mesh_c.AverageVectorFunction_perSubdomain(uv, 0); cout << "Average cup temperature: " << average_cup_temperature << " after " << time_count << " seconds. " << endl; + cout << "% of mug elements reached temperature 67ยบ: " << percentage_temp_reached << endl; time_count += dt; } auto t4 = system_clock::now(); // stop timer diff --git a/mgrid_2/visualize_results.m b/mgrid_2/visualize_results.m index 0f3546b..9295c7a 100644 --- a/mgrid_2/visualize_results.m +++ b/mgrid_2/visualize_results.m @@ -23,6 +23,7 @@ h = patch('Faces', ia, 'Vertices', xc, 'FaceVertexCData', v, ... axis equal tight; colorbar; colormap(jet); +caxis([18 85]); title('Heat distribution'); waitfor(h) % wait for closing the figure \ No newline at end of file