📄 ex18.php
字号:
</pre></div><div class = "comment">Print out status at each adaptive step.</div><div class ="fragment"><pre> Real global_error = error.l2_norm(); std::cout << "adaptive step " << a_step << ": "; if (global_tolerance != 0.) std::cout << "global_error = " << global_error << " with "; std::cout << mesh.n_active_elem() << " active elements and " << equation_systems.n_active_dofs() << " active dofs." << std::endl; if (global_tolerance != 0.) std::cout << "worst element error = " << error.maximum() << ", mean = " << error.mean() << std::endl; if (global_tolerance != 0.) {</pre></div><div class = "comment">If we've reached our desired tolerance, wedon't need any more adaptive steps</div><div class ="fragment"><pre> if (global_error < global_tolerance) break; mesh_refinement.flag_elements_by_error_tolerance(error); } else {</pre></div><div class = "comment">If flag_elements_by_nelem_target returns true, thisshould be our last adaptive step.</div><div class ="fragment"><pre> if (mesh_refinement.flag_elements_by_nelem_target(error)) { mesh_refinement.refine_and_coarsen_elements(); equation_systems.reinit(); a_step = max_adaptivesteps; break; } } </pre></div><div class = "comment">Carry out the adaptive mesh refinement/coarsening</div><div class ="fragment"><pre> mesh_refinement.refine_and_coarsen_elements(); equation_systems.reinit(); }</pre></div><div class = "comment">Do one last solve if necessary</div><div class ="fragment"><pre> if (a_step == max_adaptivesteps) { system.solve(); } </pre></div><div class = "comment">Advance to the next timestep in a transient problem</div><div class ="fragment"><pre> system.time_solver->advance_timestep(); </pre></div><div class = "comment">Write out this timestep if we're requested to</div><div class ="fragment"><pre> if ((t_step+1)%write_interval == 0) { OStringStream file_name; </pre></div><div class = "comment">We write the file name in the gmv auto-read format.</div><div class ="fragment"><pre> file_name << "out.gmv."; OSSRealzeroright(file_name,3,0, t_step + 1); GMVIO(mesh).write_equation_systems (file_name.str(), equation_systems); } } } #endif // #ifndef ENABLE_AMR </pre></div><div class = "comment">All done. </div><div class ="fragment"><pre> return libMesh::close (); }</pre></div><a name="nocomments"></a> <br><br><br> <h1> The program without comments: </h1> <pre> #include <B><FONT COLOR="#BC8F8F">"equation_systems.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"error_vector.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"getpot.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"gmv_io.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"kelly_error_estimator.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"mesh.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"mesh_generation.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"mesh_refinement.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"uniform_refinement_estimator.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"o_string_stream.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"naviersystem.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"diff_solver.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"euler_solver.h"</FONT></B> #include <B><FONT COLOR="#BC8F8F">"steady_solver.h"</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> main (<B><FONT COLOR="#228B22">int</FONT></B> argc, <B><FONT COLOR="#228B22">char</FONT></B>** argv) { <B><FONT COLOR="#5F9EA0">libMesh</FONT></B>::init (argc, argv); #ifndef ENABLE_AMR <B><FONT COLOR="#5F9EA0">std</FONT></B>::cerr << <B><FONT COLOR="#BC8F8F">"ERROR: This example requires libMesh to be\n"</FONT></B> << <B><FONT COLOR="#BC8F8F">"compiled with AMR support!"</FONT></B> << std::endl; <B><FONT COLOR="#A020F0">return</FONT></B> 0; #<B><FONT COLOR="#A020F0">else</FONT></B> { GetPot infile(<B><FONT COLOR="#BC8F8F">"ex18.in"</FONT></B>); <B><FONT COLOR="#228B22">const</FONT></B> Real global_tolerance = infile(<B><FONT COLOR="#BC8F8F">"global_tolerance"</FONT></B>, 0.); <B><FONT COLOR="#228B22">const</FONT></B> <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> nelem_target = infile(<B><FONT COLOR="#BC8F8F">"n_elements"</FONT></B>, 400); <B><FONT COLOR="#228B22">const</FONT></B> <B><FONT COLOR="#228B22">bool</FONT></B> transient = infile(<B><FONT COLOR="#BC8F8F">"transient"</FONT></B>, true); <B><FONT COLOR="#228B22">const</FONT></B> Real deltat = infile(<B><FONT COLOR="#BC8F8F">"deltat"</FONT></B>, 0.005); <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> n_timesteps = infile(<B><FONT COLOR="#BC8F8F">"n_timesteps"</FONT></B>, 20); <B><FONT COLOR="#228B22">const</FONT></B> <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> write_interval = infile(<B><FONT COLOR="#BC8F8F">"write_interval"</FONT></B>, 5); <B><FONT COLOR="#228B22">const</FONT></B> <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> coarsegridsize = infile(<B><FONT COLOR="#BC8F8F">"coarsegridsize"</FONT></B>, 1); <B><FONT COLOR="#228B22">const</FONT></B> <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> coarserefinements = infile(<B><FONT COLOR="#BC8F8F">"coarserefinements"</FONT></B>, 0); <B><FONT COLOR="#228B22">const</FONT></B> <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> max_adaptivesteps = infile(<B><FONT COLOR="#BC8F8F">"max_adaptivesteps"</FONT></B>, 10); <B><FONT COLOR="#228B22">const</FONT></B> <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> dim = infile(<B><FONT COLOR="#BC8F8F">"dimension"</FONT></B>, 2); assert (dim == 2 || dim == 3); Mesh mesh (dim); MeshRefinement mesh_refinement(mesh); mesh_refinement.coarsen_by_parents() = true; mesh_refinement.absolute_global_tolerance() = global_tolerance; mesh_refinement.nelem_target() = nelem_target; mesh_refinement.refine_fraction() = 0.3; mesh_refinement.coarsen_fraction() = 0.3; mesh_refinement.coarsen_threshold() = 0.1; <B><FONT COLOR="#A020F0">if</FONT></B> (dim == 2) <B><FONT COLOR="#5F9EA0">MeshTools</FONT></B>::Generation::build_square (mesh, coarsegridsize, coarsegridsize, 0., 1., 0., 1., QUAD9); <B><FONT COLOR="#A020F0">else</FONT></B> <B><FONT COLOR="#A020F0">if</FONT></B> (dim == 3) <B><FONT COLOR="#5F9EA0">MeshTools</FONT></B>::Generation::build_cube (mesh, coarsegridsize, coarsegridsize, coarsegridsize, 0., 1., 0., 1., 0., 1., HEX27); mesh_refinement.uniformly_refine(coarserefinements); mesh.print_info(); EquationSystems equation_systems (mesh); NavierSystem & system = equation_systems.add_system<NavierSystem> (<B><FONT COLOR="#BC8F8F">"Navier-Stokes"</FONT></B>); <B><FONT COLOR="#A020F0">if</FONT></B> (transient) system.time_solver = AutoPtr<TimeSolver>(<B><FONT COLOR="#A020F0">new</FONT></B> EulerSolver(system)); <B><FONT COLOR="#A020F0">else</FONT></B> { system.time_solver = AutoPtr<TimeSolver>(<B><FONT COLOR="#A020F0">new</FONT></B> SteadySolver(system)); assert(n_timesteps == 1); } equation_systems.init (); system.deltat = deltat; DiffSolver &solver = *(system.time_solver->diff_solver().get()); solver.quiet = infile(<B><FONT COLOR="#BC8F8F">"solver_quiet"</FONT></B>, true); solver.max_nonlinear_iterations = infile(<B><FONT COLOR="#BC8F8F">"max_nonlinear_iterations"</FONT></B>, 15); solver.relative_step_tolerance = infile(<B><FONT COLOR="#BC8F8F">"relative_step_tolerance"</FONT></B>, 1.e-3); solver.relative_residual_tolerance = infile(<B><FONT COLOR="#BC8F8F">"relative_residual_tolerance"</FONT></B>, 0.0); solver.max_linear_iterations = infile(<B><FONT COLOR="#BC8F8F">"max_linear_iterations"</FONT></B>, 50000); solver.initial_linear_tolerance = infile(<B><FONT COLOR="#BC8F8F">"initial_linear_tolerance"</FONT></B>, 1.e-3); equation_systems.print_info(); <B><FONT COLOR="#A020F0">for</FONT></B> (<B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> t_step=0; t_step != n_timesteps; ++t_step) { <B><FONT COLOR="#5F9EA0">std</FONT></B>::cout << <B><FONT COLOR="#BC8F8F">" Solving time step "</FONT></B> << t_step << <B><FONT COLOR="#BC8F8F">", time = "</FONT></B> << system.time << std::endl; <B><FONT COLOR="#228B22">unsigned</FONT></B> <B><FONT COLOR="#228B22">int</FONT></B> a_step = 0; <B><FONT COLOR="#A020F0">for</FONT></B> (; a_step != max_adaptivesteps; ++a_step) { system.solve(); ErrorVector error; AutoPtr<ErrorEstimator> error_estimator; <B><FONT COLOR="#A020F0">if</FONT></B> (global_tolerance != 0.) { assert (nelem_target == 0); UniformRefinementEstimator *u = <B><FONT COLOR="#A020F0">new</FONT></B> UniformRefinementEstimator; u->sobolev_order() = 0; error_estimator.reset(u); } <B><FONT COLOR="#A020F0">else</FONT></B> { assert (nelem_target > 0); error_estimator.reset(<B><FONT COLOR="#A020F0">new</FONT></B> KellyErrorEstimator); } error_estimator->component_scale.push_back(1.0); <I><FONT COLOR="#B22222">// u</FONT></I> error_estimator->component_scale.push_back(1.0); <I><FONT COLOR="#B22222">// v</FONT></I> <B><FONT COLOR="#A020F0">if</FONT></B> (dim == 3) error_estimator->component_scale.push_back(1.0); <I><FONT COLOR="#B22222">// w</FONT></I> error_estimator->component_scale.push_back(0.0); <I><FONT COLOR="#B22222">// p</FONT></I> error_estimator->estimate_error(system, error); Real global_error = error.l2_norm(); <B><FONT COLOR="#5F9EA0">std</FONT></B>::cout << <B><FONT COLOR="#BC8F8F">"adaptive step "</FONT></B> << a_step << <B><FONT COLOR="#BC8F8F">": "</FONT></B>; <B><FONT COLOR="#A020F0">if</FONT></B> (global_tolerance != 0.) <B><FONT COLOR="#5F9EA0">std</FONT></B>::cout << <B><FONT COLOR="#BC8F8F">"global_error = "</FONT></B> << global_error << <B><FONT COLOR="#BC8F8F">" with "</FONT></B>; <B><FONT COLOR="#5F9EA0">std</FONT></B>::cout << mesh.n_active_elem() << <B><FONT COLOR="#BC8F8F">" active elements and "</FONT></B> << equation_systems.n_active_dofs() << <B><FONT COLOR="#BC8F8F">" active dofs."</FONT></B> << std::endl; <B><FONT COLOR="#A020F0">if</FONT></B> (global_tolerance != 0.) <B><FONT COLOR="#5F9EA0">std</FONT></B>::cout << <B><FONT COLOR="#BC8F8F">"worst element error = "</FONT></B> << error.maximum() << <B><FONT COLOR="#BC8F8F">", mean = "</FONT></B> << error.mean() << std::endl; <B><FONT COLOR="#A020F0">if</FONT></B> (global_tolerance != 0.) { <B><FONT COLOR="#A020F0">if</FONT></B> (global_error < global_tolerance) <B><FONT COLOR="#A020F0">break</FONT></B>; mesh_refinement.flag_elements_by_error_tolerance(error); } <B><FONT COLOR="#A020F0">else</FONT></B> { <B><FONT COLOR="#A020F0">if</FONT></B> (mesh_refinement.flag_elements_by_nelem_target(error)) { mesh_refinement.refine_and_coarsen_elements(); equation_systems.reinit(); a_step = max_adaptivesteps; <B><FONT COLOR="#A020F0">break</FONT></B>; } } mesh_refinement.refine_and_coarsen_elements(); equation_systems.reinit(); } <B><FONT COLOR="#A020F0">if</FONT></B> (a_step == max_adaptivesteps) { system.solve(); } system.time_solver->advance_timestep(); <B><FONT COLOR="#A020F0">if</FONT></B> ((t_step+1)%write_interval == 0) { OStringStream file_name; file_name << <B><FONT COLOR="#BC8F8F">"out.gmv."</FONT></B>; OSSRealzeroright(file_name,3,0, t_step + 1); GMVIO(mesh).write_equation_systems (file_name.str(), equation_systems); } } } #endif <I><FONT COLOR="#B22222">// #ifndef ENABLE_AMR</FONT></I> <B><FONT COLOR="#A020F0">return</FONT></B> libMesh::close (); }</pre> <a name="output"></a> <br><br><br> <h1> The console output of the program: </h1> <pre>**************************************************************** Running Example ./ex18-devel*************************************************************** Mesh Information: mesh_dimension()=2 spatial_dimension()=3 n_nodes()=1681 n_elem()=400 n_local_elem()=400 n_active_elem()=400 n_subdomains()=1 n_processors()=1 processor_id()=0*** Warning, This code is untested, experimental, or likely to see future API changes: src/solvers/diff_system.C, line 29, compiled Jun 1 2007 at 14:30:46 *** EquationSystems n_systems()=1 System "Navier-Stokes" Type "Implicit" Variables="u" "v" "p" Finite Element Types="LAGRANGE" "LAGRANGE" "LAGRANGE" Approximation Orders="SECOND" "SECOND" "FIRST" n_dofs()=3803 n_local_dofs()=3803 n_constrained_dofs()=0 n_vectors()=2 Solving time step 0, time = 0 Solving time step 1, time = 0.005 Solving time step 2, time = 0.01 Solving time step 3, time = 0.015 Solving time step 4, time = 0.02 Solving time step 5, time = 0.025 Solving time step 6, time = 0.03 Solving time step 7, time = 0.035 Solving time step 8, time = 0.04 Solving time step 9, time = 0.045 Solving time step 10, time = 0.05 Solving time step 11, time = 0.055 Solving time step 12, time = 0.06 Solving time step 13, time = 0.065 Solving time step 14, time = 0.07 **************************************************************** Done Running Example ./ex18-devel***************************************************************</pre></div><?php make_footer() ?></body></html><?php if (0) { ?>\#Local Variables:\#mode: html\#End:<?php } ?>
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -