Thank you for looking at the code, I thought that the filter would do something magic (I was naive), but of course it cannot defied logical rules.
I ended up refactoring most of the code, looping first over the active cells, and now it is going much better. On Tue, May 28, 2024, 7:13 PM Wolfgang Bangerth <[email protected]> wrote: > On 5/27/24 00:28, Felipe Ponce wrote: > > > > Probably, the problem is not the filtered iterator, but maybe I should > > re-think the algorithm, to avoid a costly search for active descendants. > > Yes. The three methods you show in your code are doing very different > things: > > // Method I: Iterate using indices. > for (index_t i = 0; i < coarse_cell.n_children(); ++i) > { > const index_t id_child = coarse_cell.child(i)->index(); > } > > Here, for each of the 384 coarse mesh cells in the outer loop (not shown) > you > are only looping over the 4 children. So you are touching 384*4 = 1536 > cells. > > > // Method II: Iterate using filter. > ChildIt cell(ActiveChild(coarse_cell.id())); > ChildIt endc(ActiveChild(coarse_cell.id()), dof_handler.end()); > cell.set_to_next_positive(dof_handler.begin_active()); > for (; cell != endc; ++cell) > { > cell->get_dof_indices(global_dofs); > } > > Here, for each of the 384 coarse mesh cells, you are looping over *all* > 1536 > fine mesh cells. So you are touching 384*1536 = 589,824 cells in either > the > set_to_next_positive() call or in the ++cell calls where you need to test > whether a cell is a descendant of the current coarse cell. > > > // Method II: Iterate without using filter. > for (const auto &cell : dof_handler.active_cell_iterators()) > { > if (coarse_cell.id().is_ancestor_of(cell->id())) > { > cell->get_dof_indices(global_dofs); > } > } > > And here you do the same, just by hand. > > So you are doing fundamentally different things. I'm not surprised you get > fundamentally different run times :-) > > Best > W. > > -- > ------------------------------------------------------------------------ > Wolfgang Bangerth email: [email protected] > www: http://www.math.colostate.edu/~bangerth/ > > > -- > 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 [email protected]. > To view this discussion on the web visit > https://groups.google.com/d/msgid/dealii/c501a2b0-69c9-430d-87d0-eef1d765ab35%40colostate.edu > . > -- 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 [email protected]. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/CAEjcB7VO2K-AdNmiy6fJWbHsLMCfFU%2BU71AZXyz94k-7%2BXL9FQ%40mail.gmail.com.
