Dear Prof. Wolfgang, Thank you very much for your quick answer! And I have tried what you say, and I have some questions still: 1.
> In your case, a > SparsityPattern is still used by some other object but is being destroyed. > You'll have to find out what other object that is -- likely a SparseMatrix > object. You'll first have to break that dependence before you destroy the > sparsity pattern. > I had read : https://groups.google.com/forum/#!searchin/dealii/Object$20of$20class$20N6dealii15SparsityPatternE$20is$20still$20used$20by$201$20other$20objects%7Csort:date/dealii/8Sqkbo1n5Jg/O6R2Vpya4hwJ And I think what you means is the same as it above beforehand, but maybe not now. The dsp (BlockDynamicSparsityPattern) will be destroyed and cannot exist in all lifetime of BlockSparseMatrix, so we define : BlockSparsityPattern sparsity_pattern at the top, (this is what step-22 does) and let it copy dsp.... But now, we know that BlockSparsityPattern sparsity_pattern also will be destroyed and cannot exist in all lifetime of BlockSparseMatrix. So I don't know how to correct it. And I cannot implement what you said "You'll first have to break that dependence before you destroy the sparsity pattern". 2. As your suggestions, I change my code : template <int dim> void StokesProblem<dim>::refine_mesh (const unsigned int max_grid_level, const unsigned int min_grid_level) { Vector<float> estimated_error_per_cell (triangulation_active.n_active_cells()); FEValuesExtractors::Scalar phase(0); KellyErrorEstimator<dim>::estimate (dof_handler, QGauss<dim-1>(degree+1), typename FunctionMap<dim>::type(), old_solution, estimated_error_per_cell, fe.component_mask(phase)); GridRefinement::refine_and_coarsen_fixed_number (triangulation_active, estimated_error_per_cell, 0.2,0.1); std::cout << "triangulation.n_levels" << triangulation_active.n_levels() << std::endl; if (triangulation_active.n_levels() > max_grid_level) for (typename Triangulation<dim>::active_cell_iterator cell = triangulation.begin_active(max_grid_level); cell != triangulation.end(); ++cell) cell->clear_refine_flag (); std::cout << "clear-refine" << triangulation_active.n_levels() << std::endl; if (triangulation_active.n_levels() < min_grid_level) for (typename Triangulation<dim>::active_cell_iterator cell = triangulation.end_active(min_grid_level); cell != triangulation.end(); ++cell) cell->clear_coarsen_flag (); std::cout << "clear-coarsen" << triangulation_active.n_levels() << std::endl; triangulation_active.execute_coarsening_and_refinement (); std::cout << "c-and-r" << triangulation_active.n_levels() << std::endl; } Maybe I need to chage all "triangulation" in my originial StokesProblem<dim>::refine_mesh (const unsigned int max_grid_level, const unsigned int min_grid_level) into "triangulation_active", this function can run(but i don't know if it is right), and for implement periodic boundary condition, I use the " triangulation" and "triangulation_active" : for (unsigned int cycle=0; cycle<n_cycles; ++cycle) { if (cycle == 0) { GridGenerator::hyper_cube (triangulation_active, 0, 4.5); triangulation_active.refine_global (8); std::cout << "first-triangulation.n_levels" << triangulation_active.n_levels() << std::endl; GridGenerator::flatten_triangulation (triangulation_active, triangulation); std::cout << "second-triangulation.n_levels" << triangulation.n_levels() << std::endl; for (typename Triangulation<dim>::active_cell_iterator cell = triangulation.begin_active(); cell != triangulation.end(); ++cell) { for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f) { if (cell->face(f)->at_boundary()) { if (std::fabs(cell->face(f)->center()(0) - (4.5)) < 1e-12) cell->face(f)->set_boundary_id (1); else if(std::fabs(cell->face(f)->center()(1) - (0)) < 1e-12) cell->face(f)->set_boundary_id (2); else if(std::fabs(cell->face(f)->center()(1) - (4.5)) < 1e-12) cell->face(f)->set_boundary_id (3); else if(std::fabs(cell->face(f)->center()(0) - (0)) < 1e-12) { cell->face(f)->set_boundary_id (0); } } } } } I don't know what I should use in whole code is "triangulation" or " triangulation_active" ? And before I add the adaptive mesh, the original code can run, and I used "triangulation_active" only the definition of mesh at first, and "triangulation" is used after all these. But now, I am confused that maybe I need to use "triangulation_active" in refine_mesh() function, if it works, I donnot know the reason. Hope your answer and thank you very much! Best, chucui -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. For more options, visit https://groups.google.com/d/optout.