sandialabs
Limiting the number of elements
Good afternoon, I am an engineering student working on a project on mesh adaptation and I am using Omega_h for the first time. I have to deal with a 3d mesh of tetrahedrons that I wish to adapt in order to underline the presence of a ball-shaped figure. I have constructed the mesh using the program Fenics and computed the metric manually before using Omega_h. The whole program is structured in 2 nested while loops plus the one of adaptation in Omega_h. The adaptation seems to be going in the right direction as it is possible to see in the picture below, however I am having issues dealing with the number of elements of the mesh. In particular the initial (isotropic) mesh was of 39 546 tetrahedrons, after the first cycle of "while (Omega_h::approach_metric...", which lasts 6 iterations, the number of elements increases to 640 599 tetrahedrons , which is I believe acceptable. The problem is that the second time i call the adaptation procedure (again it lasts 6 iterations) the number of elements increases to 1.17921e+07 elements, which is way too much. To solve the problem i first tried to add the lines: metric_input.should_limit_element_count = true; metric_input.max_element_count = 2000000 But this resulted in the code blocking (without errors) at the line " Omega_h::generate_target_metric_tag(&mesh_osh, metric_input);". I have also tried to change some parameters but without success, both in the metric construction part of the code and in the Omega_h part.
For completeness this is my code for the function that adapts the mesh with Omega_h:
`void Adaptation(Omega_h::Library & lib_osh, std::shared_ptr<Mesh> Th,
const std::vector
Omega_h::Mesh mesh_osh(&lib_osh); //Initialize the Omega_h metric Omega_h::from_dolfin(&mesh_osh, Th); //Transfom the dolfin mesh in Omega_h mesh Omega_h::classify_by_angles(&mesh_osh, Omega_h::PI / 4.0);//Labelize every element of the mesh Omega_h::Write<Omega_h::Real> vs(6mesh_osh.nverts()); //Initialize writing object (dummy variable that will contain our metric) Omega_h::MetricInput metric_input; //Initialize Input Metric (Will be used for adaptation) Omega_h::add_implied_metric_tag(&mesh_osh); //Add metric tags to the Omega_h mesh
for( int i=0; i<mesh_osh.nverts();i++) { Omega_h::Vector<6> v = {m11[i], m22[i], m33[i], m12[i], m23[i], m13[i]}; Omega_h::set_vector(vs, i, v); //Set the long vector containing piecewise metric vectors }
Omega_h::Read<Omega_h::Real> vc=vs; //Transform the write object into a read object mesh_osh.add_tag(Omega_h::VERT, "new", Omega_h::symm_ncomps(mesh_osh.dim()), vc); //Add the new metric tag to the mesh auto source = Omega_h::MetricSource(OMEGA_H_GIVEN, 1.0, "new"); //Intialize the Object containing new metric metric_input.sources.push_back(source); //Input metric is the new metric metric_input.should_limit_gradation = true; //Limit the gradation (see reference for this) metric_input.max_gradation_rate=5.0;
Omega_h::vtk::write_vtu("/results/beforegenerating.vtu", &mesh_osh); //Save the old mesh in a vtu (before even generating the metric) Omega_h::generate_target_metric_tag(&mesh_osh, metric_input); //Insert the tag of the new metric Omega_h::vtk::write_vtu("/results/before.vtu", &mesh_osh); //Save the old mesh but with the new Metrics tags Omega_h::AdaptOpts opts(&mesh_osh); //Initialize Options for the adaptation opts.min_quality_allowed=1e-6; //Minimum Flattening of the triangles (related to Anisotropy)
opts.max_length_allowed=3.2;
Omega_h::vtk::Writer writer("/results/adapting", &mesh_osh); int count=0; while (Omega_h::approach_metric(&mesh_osh, opts)) { std::cout<<"While di Omega_H numero: " << count << std::endl; count++; Omega_h::adapt(&mesh_osh, opts); writer.write(); } Omega_h::vtk::write_vtu("/results/after.vtu", &mesh_osh); Omega_h::to_dolfin(*Th, &mesh_osh); Th->order(); `
I was wondering if there is some clear mistake I may be doing or if since my code already contains a while loop that checks the error and the variation of the number of elements, I could avoid the "while (Omega_h::approach_metric..." part. Moreover I would like to understand better the parameters in "AdaptOpts::AdaptOpts(Int dim)" which is a function of Omega_h_adapt.cpp. To make the code work I had to in fact manually add the line "opts.max_length_allowed=3.2".
Thank you very much in advanced to anyone who will help me. Sharon Dacatra
