Iterating over two DoFhandlers for a thermo-elastic problem

[subramanya gautam]

I am trying to implement a staggered thermal-elastic solver using the matrix-free multigrid solver approach from step-37. I got an elasticity solver and a thermal solver working and now I am trying to associate the two. I have a shared triangulation, with two separate dof_handlers for the  thermal and stress problems. The question I have is when I do the assembly loop, does the cell refer to the same cell for both the thermal and  elastic problems? 

So from step-37,  can I just call phi_temp.reinit(cell), and phi_stress.reinit(cell), and have cell refer to the same physical cell in the triangulation? 

for (unsigned int cell = 0;

  cell < system_matrix.get_matrix_free()->n_cell_batches();

  ++cell)

  {

  phi.reinit(cell);

  for (const unsigned int q : phi.quadrature_point_indices())

  phi.submit_value(make_vectorized_array<double>(1.0), q);

  phi.integrate(EvaluationFlags::values);

  phi.distribute_local_to_global(system_rhs);

  }

Thanks!

Subramanya.

[Wolfgang Bangerth]

On 1/3/25 09:14, subramanya gautam wrote:
> I am trying to implement a staggered thermal-elastic solver using the
> matrix-free multigrid solver approach from step-37. I got an elasticity solver
> and a thermal solver working and now I am trying to associate the two. I have
> a shared triangulation, with two separate dof_handlers for the  thermal and
> stress problems. The question I have is when I do the assembly loop, does the
> cell refer to the same cell for both the thermal and  elastic problems?
> So from step-37,  can I just call phi_temp.reinit(cell), and
> phi_stress.reinit(cell), and have cell refer to the same physical cell in the
> triangulation?

Yes, the order of cells of DoFHandlers follows the order of the triangulation.
So two DoFHandlers using the same triangulation have the same order of cells.
step-31, for example, uses this fact.

Best
Wolfgang

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

[Subramanya G]

Thank you so much! 

Subramanya.  

सुब्रह्मण्य .

--
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+un...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/dealii/f40ad605-daf4-4347-b4ff-2b589c5b04cb%40colostate.edu.

[Martin Kronbichler]

Dear Subramanya,

The actual cells associated with the integer indices 'cell' (that are batches in case of SIMD vectorization) are not necessarily the same if you create two separate MatrixFree objects and initialize each with the respective DoFHandler for the two fields. If you intend to use coupled problems like the thermal-stress problem you describe, you should use the MatrixFree::reinit function that takes both DoFHandler objects as the arguments to an std::vector<DoFHandler*> and respectively for the affine constraints. Then, you can initialize FEEvaluation objects using the 'dof_no' argument in the constructor: https://www.dealii.org/developer/doxygen/deal.II/classFEEvaluation.html#ade40e452949cbefd4f7e6d5580af81c8 - it is the second argument next to MatrixFree.

You can also try to run with separate MatrixFree objects, but then you need to manually link the two objects and one FEEvaluation needs to be initialized with an array of cells (that you constructed manually). I have done that in some applications, and there is at least one code inside deal.II that does this: https://github.com/dealii/dealii/blob/8f97fd0dca1eb878bc7ae9922ac89430cab800db/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h#L4041-L4068 - but it is much less obvious, so I really recommend the first variant.

I hope this helps!

Best regards,
Martin

Am 03.01.25 um 17:14 schrieb subramanya gautam:

--
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+un...@googlegroups.com.

To view this discussion visit https://groups.google.com/d/msgid/dealii/d8f59b24-2a30-4ade-83e4-75bc72112afen%40googlegroups.com.

[Subramanya G]

Thanks Martin for all the help! 

If I understand correctly, the first approach would assume that both the problems are solved in a fully coupled manner, with all fields updated simultaneously at each solve. 

For a staggered solution , where I solve the temperature first, and then the stress - I'd either have to manually compute the cell matching between the two objects, before assembling the vector, or abandon the SIMD vectorization, and compute the vector the "old-fashioned" way just for the coupling where information from the two vectors needs to be accessed simultaneously.  

Subramanya.  

सुब्रह्मण्य .

To view this discussion visit https://groups.google.com/d/msgid/dealii/71ea7f4c-19c8-43f6-9037-775aa504b470%40rub.de.

[Martin Kronbichler]

Dear Subramanya,

The first option does not tie you to solving fully coupled problems. In fact, you can do any kind of staggered/coupled approach as you want, because MatrixFree is conceptually just the data container collecting various things, whereas you then define one or more operators on top of it, like the LaplaceOperator in step-37. So when you go with the first approach (what I recommend), you first set up a single MatrixFree object with both DoFHandler objects (and all other relevant information, you can also give different quadrature formulas, for example). Then you create two classes, which are both using the same MatrixFree object (or via the same shared_ptr), where the vmult() function in one class uses "FEEvaluation<dim,....> eval_stress(matrix_free, 0)" as the evaluator object and the vmult() function in the other class uses "FEEvaluation<dim,....> eval_temp(matrix_free, 1)". And then you have a coupling routine (organized in any way you like, but again with access to the main MatrixFree object) that has both these FEEvaluation objects in its evaluation function and does whatever you want to do, e.g., extract the values from one field (via read_dof_values()/evaluate()) and consuming them as "coefficient" in the other field (via integrate()/distribute_local_to_global()) on the other field, which happens at the quadrature points.

As long as they are in the same MatrixFree object, everything is SIMD vectorized and assembled in an optimized way.

Best regards,
Martin

Am 05.01.25 um 04:30 schrieb Subramanya G:

To view this discussion visit https://groups.google.com/d/msgid/dealii/CAHCJNxN82c0u_L8VU31gBS0tK8WxpxkJXmjwB%2Bav5P2jx%2B4-qw%40mail.gmail.com.