Hi, Dr. W. Bangerth and everyone,
Thank you for your help, Dr. W. Bangerth!
I have figured the solution out. I summarize it as follows in case that it
can help some one else as well.
*Problem summary:* In topology optimization, each cell has a scalar
variable (pesudo density). So the density will form a vector, but the
vector dimensions are different from those of the vector, say
LA::MPI::Vector locally_relevant_solution (refer to Step-40). Because it
only has one related value per cell. How to build up such "cell-based"
vector?
*Solution steps:* (inspired by the suggestion from Dr. W. Bangerth)
1. Build two fem field objects.
FESystem<dim> fe; // for describing the vector-valued field (eg.
displacement)
FE_DGQ<dim> fe_dg; // for describing the cell field (eg. pesudo density in
topology optimization)
And initialize them in the initialization list of the class,
fe(FE_Q<dim>(1), dim),
fe_dg(0),
And initialize the cell-based vector (the pesudo density in topology
optimization) as follows,
opt->locally_owned_cells= opt->dof_handler_dg.locally_owned_dofs();
// it is named with "cells" because I want to distinguish it from
"opt->locally_owned_dofs" which is for the first FE field (fe(FE_Q<dim>(1),
dim)).
DoFTools::extract_locally_relevant_dofs(opt->dof_handler_dg,
opt->locally_relevant_cells);
opt->x.reinit(opt->locally_owned_cells,
opt->locally_relevant_cells,
opt->mpi_communicator);
2. I need to access both FE fields at the same time.
DoFHandler<dim> dof_handler;
DoFHandler<dim> dof_handler_dg;
If I use the default loop,
for (const auto &cell : opt->dof_handler.active_cell_iterators())
{
...
}
how can I access to the second FE field, which is dof_handler_dg?
My attempting solution is to set up multiple iterators, just as in the
tutorial:
https://www.dealii.org/current/doxygen/deal.II/group__Iterators.html,
typename Triangulation<dim>::active_cell_iterator ti =
opt->triangulation.begin_active();
typename DoFHandler<dim>::active_cell_iterator di1 =
opt->dof_handler.begin_active();
typename DoFHandler<dim>::active_cell_iterator di2 =
opt->dof_handler_dg.begin_active();
Therefore, I set up the loop as follows,
while (cell != opt->triangulation.end())
{
// do someting
++ti;
++di1;
++di2;
}
3. Now I can access to the cell-based vector x successfully as follows,
di2->get_dof_indices(local_dof_indices_dg);
unsigned int local_dof_indices_dg_int = local_dof_indices_dg.at(0);
opt->x(local_dof_indices_dg_int);
It was very time-consuming to execute because I put the vector access, say
opt->x(local_dof_indices_dg_int), in the most inner loop. I don't need to
access it for each quadrature point nor each degree of freedom. Because for
each cell, the opt->x is fixed. So I put opt->x in the outside of those
loops, and the computing time is normal now.
In sum, the parallel distributed density vector has been set up.
Thank you very much for Dr. W. Bangerth's help!
Best regards,
Z. B. Zhang
On Sunday, October 20, 2019 at 9:26:59 PM UTC-4, Wolfgang Bangerth wrote:
>
>
> > > Yes. Though it seems confusing to me that you call these variables
> > > locally_*_cells, because they really are index sets for degrees of
> freedom,
> > > not cells.
> >
> > OK, I misunderstood it. I thought the numbering of the dof is same as
> that of
> > the cell since the degrees of freedom is 1 in FE_DGQ(0).
> > They are still different, right?
>
> They are *conceptually* different. Whether or not they are different *in
> practice* is a different question, but you should never *assume* that they
> are
> the same.
>
>
> > > Yes, this does not work. That's because you confuse the index of the
> cell and
> > > the index of the DoF that lives on it. You want the latter.
> >
> > Again now I understand they are different.
> >
> > However, for my project, I need to access both
> >
> > DoFHandler<dim> dof_handler;
> > DoFHandler<dim> dof_handler_dg;
> >
> > at the same time.
> >
> > If I use the following loop,
> >
> > for (const auto &cell : opt->dof_handler.active_cell_iterators())
> > {
> > ...
> > }
> >
> > how can I access to dof_handler_dg?
>
> You can have two iterators into the two DoFHandlers that run in lock step
> just
> like you do:
>
>
> > My attempting solution is to set up multiple iterator, just as in the
> tutorial:
> > https://www.dealii.org/current/doxygen/deal.II/group__Iterators.html,
> >
> > typename Triangulation<dim>::active_cell_iterator ti =
> > opt->triangulation.begin_active();
> > typename DoFHandler<dim>::active_cell_iterator di1 =
> > opt->dof_handler.begin_active();
> > typename DoFHandler<dim>::active_cell_iterator di2 =
> > opt->dof_handler_dg.begin_active();
> >
> > Therefore, I set up the loop as follows,
> >
> > while (cell != opt->triangulation.end())
> > {
> > // do someting
> >
> > ++ti;
> > ++di1;
> > ++di2;
> > }
> >
> > (Note: now I can access to the cell-based vector x successfully as
> follows,
> > di2->get_dof_indices(local_dof_indices_dg);
> > unsigned int local_dof_indices_dg_int = local_dof_indices_dg.at(0);
> > opt->x(local_dof_indices_dg_int);
> > )
> >
> >
> > However, I found this method is very time-consuming compared to "for
> (const
> > auto &cell : opt->dof_handler.active_cell_iterators())", like 11s vs. 1
> s.
>
> You mean time consuming to program, or to execute? The loop itself is
> cheap;
> if it's slow, you need to figure out which part of the body is slow.
>
>
> > My question now is: is there other method to access multiple DoFHandler
> or/and
> > Triangulation?
>
> You could just do everything in one DoFHandler. That's what I usually do:
> Put
> all solution variables into the same DoFHandler.
>
> Best
> W.
>
> --
> ------------------------------------------------------------------------
> Wolfgang Bangerth email: [email protected]
> <javascript:>
> 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/8e2bc60b-4779-450e-9322-a374e7965e9b%40googlegroups.com.
