How can I visit a face exactly once when running with MPI?

63 views
Skip to first unread message

Abbas

unread,
Apr 27, 2023, 7:56:50 AM4/27/23
to deal.II User Group
I want to postporcess some data across an interface with feInterface. 
I have a fully::distributed::triangulation. 
I am doing is something that "roughly" looks like this: 

template <int dim>
void TSL_prllel<dim>::interface_gradients(Tensor<1, dim> &interface_gradients)
{
interface_gradients = 0;
QGauss<dim - 1> face_quadrature_formula(degree + 1);

FEInterfaceValues<dim> fe_iv(fe,
face_quadrature_formula,
+many_flags);

Tensor<1, dim> local_interface_gradients;

for (; cell != endc; ++cell)
if (cell->is_locally_owned())
for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell; ++face)
if (cell->face(face)->user_flag_set()) //If I am at the interface
{
// I want to visit the face only once so I clear the flag
cell->face(face)->clear_user_flag();

fe_iv.get_jump_in_function_values(ghosted_solution, gradu);
for (unsigned int point = 0; point < n_face_q_points; ++point)
{

local_interface_gradients += gradu[point]xJXW ;

}
}
interface_stress = Utilities::MPI::sum(local_interface_gradients, mpi_communicator);
}

I basically set a flag at the faces I want to postprocess at then when I visit the face the first time, I clear the flag so the face gets visited only once. This works with a single MPI but not for more. Is it because the cell->face(face)->clear_user_flag()
doesn't clear the user flag in the other distributed meshes?
I can't directly loop over faces because there is no "face->is_locally_owned()"
I know that meshworker is an option. Is it the only option?

Bruno Turcksin

unread,
Apr 27, 2023, 8:48:23 AM4/27/23
to deal.II User Group
Hello,

On Thursday, April 27, 2023 at 7:56:50 AM UTC-4 abbas.b...@gmail.com wrote:
This works with a single MPI but not for more. Is it because the cell->face(face)->clear_user_flag()
doesn't clear the user flag in the other distributed meshes?

That's correct. The flag is only cleared locally not on the other processor.
I can't directly loop over faces because there is no "face->is_locally_owned()"
I know that meshworker is an option. Is it the only option?
You can check if the cell that shares the face (use this) is a ghost cell (use this). If that's the case, you can ask for the subdomain id of the locally owned cell and the ghost cell (use this). Now you can decide that only the processor with the lower rank can "touch" the face and the other processor should just skip it.

Best,

Bruno
 

Abbas Ballout

unread,
Apr 27, 2023, 10:09:54 AM4/27/23
to deal.II User Group
This works!! Thanks. 
I am guessing something like: 

if(cell_neighbour->is_ghost()==true)
cell->face(neighbour_face_number)->clear_user_flag();

won't clear the user flag everywhere either, because that didn't work. 
Right?  

Bruno Turcksin

unread,
Apr 27, 2023, 11:18:14 AM4/27/23
to dea...@googlegroups.com
Right, it won't clear it everywhere. You would need to use MPI to clear the flag on the other processor. If you would do that you would create a huge amount of messages and your code would be extremely slow.

Best,

Bruno

--
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 a topic in the Google Groups "deal.II User Group" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/dealii/s07h9mCZhq4/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dealii+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/d7b3b433-e421-4b2c-a837-add3bd1d6f2en%40googlegroups.com.

Timo Heister

unread,
Apr 27, 2023, 5:29:34 PM4/27/23
to dea...@googlegroups.com
I would suggest that you use MeshWorker::mesh_loop() to do this (or
look at the complicated implementation of it).

On Thu, Apr 27, 2023 at 11:18 AM Bruno Turcksin
<bruno.t...@gmail.com> wrote:
>
> Right, it won't clear it everywhere. You would need to use MPI to clear the flag on the other processor. If you would do that you would create a huge amount of messages and your code would be extremely slow. Best, Bruno ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍ ‍
> ZjQcmQRYFpfptBannerStart
> This Message Is From an External Sender
> Use caution when opening links or attachments if you do not recognize the sender.
>
> ZjQcmQRYFpfptBannerEnd
> 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+un...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/CAGVt9ePNjo-9oOPy5YdfZ-RntzJDe7eu5Zdeyy6vFL8W9Jod_A%40mail.gmail.com.



--
Timo Heister
http://www.math.clemson.edu/~heister/

Abbas Ballout

unread,
Jun 7, 2023, 3:21:56 AM6/7/23
to deal.II User Group
What if I am setting the user flag with something like this: 

for (const auto &cell : triangulation.active_cell_iterators())
for (const unsigned int face : cell->face_indices())
if ((cell->face(face)->center is between x1 and x2 )
cell->face(face)->set_user_flag();

I set the flag when I create the mesh, then I loop on the faces inside a mesh worker mesh loop and do something special at a flagged face.
This has been working without any issues for sometime now. Is this error prone?  
 

Wolfgang Bangerth

unread,
Jun 7, 2023, 10:40:51 AM6/7/23
to dea...@googlegroups.com
On 6/7/23 01:21, Abbas Ballout wrote:
> I set the flag when I create the mesh, then I loop on the faces inside a mesh
> worker mesh loop and do something special at a flagged face.
> This has been working without any issues for sometime now. Is this error prone?

This works on a single process. But it will not work if you are running in
parallel because clearing the flag for a face on one process does not clear
the flag for that face on another process.

Best
W.

--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/


Abbas Ballout

unread,
Jun 7, 2023, 2:02:31 PM6/7/23
to deal.II User Group
Thanks Wolfgang, Time, and Bruno for your comments. I hate to say thanks too much so sometimes I don't. 
Reply all
Reply to author
Forward
0 new messages