Transferring gradients between submeshes for Neumann boundary conditions

[bclyons12]

In the multidomain miniapp, there is a clear example for how to transfer the value of a solution from one submesh to another on a boundary in order to impose Dirichlet boundary conditions. You make a transfer map, then perform the transfer:

auto temperature_block_to_cylinder_map = ParSubMesh::CreateTransferMap(
                                               temperature_block_gf,
                                               temperature_cylinder_gf);

then

temperature_block_to_cylinder_map.Transfer(temperature_block_gf, temperature_cylinder_gf);

What if instead I want to know the gradient on the boundary from the first submesh and use that as a Neumann boundary condition in the second submesh? Is there any example of this?

[jandrej]

There is no clear example of this as it is very problem specific. The process should go like

• Find quadrature point locations on primary side (the locations where we need the QoI)
• Find the elements on the secondary mesh where we can interpolate those coordinates
• Compute QoI on secondary elements and save result in primary "quadrature point data"
• Use this "quadrature point data" on the primary mesh to compute your FE function

See https://github.com/mfem/mfem/blob/navier-urans/miniapps/navier/tests/test_normal_stress.cpp.

This works fine in serial but needs some kind of communication scheme in parallel. This is currently work in progress to be a feature with gslib

Hope this helps.

[bclyons12]

Thanks @jandrej. I see it's trickier than I thought. The typical Neumann condition requires that the appropriate derivative is a constant along the boundary. We need to do something like this more generic method you've posted so that the derivative is the correct value at each point.

[colmenaresj]

Hi @jandrej , thanks for the clarification. You mentioned that the code works well in serial, but needs some work for it to work in parallel. Looking at the code, I see the code uses ParMesh and the ContactQoiEvaluator seems to be using GSLIB (through the FindPointsGSLIB object). Could you please explain what pieces are missing for it to work in parallel?

[jandrej]

If the mesh is actually distributed, the element numbers don't match on each MPI rank. That part of the algorithm is basically missing. gslib knows which element it found on what MPI rank but the gather/scatter requires a custom data package to be communicated. We're currently trying to get gslib to accommodate that.

[colmenaresj]

thank you for clarifying

[kmittal2]

@colmenaresj @bclyons12 I have setup #4270 that implements the functionality you are looking for. Please test it out for use case.

[colmenaresj]

Hi @kmittal2 , I have some time now to give this a try. However, for us it is important that we can run on GPUs. Does GSLIB work on GPUs? Or of it has to run on CPU, will it add a big overhead?
Thanks.

[kmittal2]

@colmenaresj FindPointsGSLIB for GPUs is currently under development, so running on CPUs is the only option. The overhead depends on your problem size.