Time dependent Step 46, implementing non-homogeneous kinetic coupling

102 views
Skip to first unread message

Sclah

unread,
Dec 17, 2024, 3:27:49 PM12/17/24
to deal.II User Group
Hello everyone,

I am using step-46 as starting point to implement a time dependent fsi problem using block vectors and block matrices.
I successfully written the time iteration and changed everything to block structures but I am struggling to find a method to impose the fluid velocity equal to the displacement velocity at the interface (i.e. setting a fsi kinetic coupling).
My approach is to use AffineConstraints::add_entry and AffineConstraints::set_inhomogeneity to write something like v=(u - u_old) / time_step, where u_old is the previous time step solution.
The problem is to find the right column index j required by these class methods. Is there a way to find it?
Thank you in advance
 

Wolfgang Bangerth

unread,
Dec 17, 2024, 10:27:19 PM12/17/24
to dea...@googlegroups.com
Where do you need them, specifically? Can you use
cell->face(f)->get_dof_indices(...)
for example to query the DoF indices on the face of a cell, where you're only
doing this on faces that separate the two kinds of materials?

Best
W.

Message has been deleted

Sclah

unread,
Dec 18, 2024, 6:00:05 PM12/18/24
to deal.II User Group
I need them just on the faces that separate the materials (fluid and solid). The piece of code is the following, on an interface cell face I have:

cell->face(face_number)->get_dof_indices(local_face_dof_indices_fluid, fluid_fe_index); // select the fluid dofs sets from hp-dofs
cell->face(face_number)->get_dof_indices(local_face_dof_indices_solid, solid_fe_index); // select the solid dofs sets from hp-dofs

// impose fluid_velocity = solid_velocity (where solid_velocity = (solid_displacement - solid_displacement_previous_time_step) / time_step)

for (unsigned int i = 0; i < local_face_dof_indices_fluid.size(); ++i)
{
if (fe_fluid->face_system_to_component_index(i).first < dim) // select just fluid velocity (the solution is stored as (fluid_velocity + fluid_pressure + solid_displacement) with (dim + 1 + dim) components)
{
auto column_index = /*todo*/ ; // need to extract the right index to select the displacement value at the interface
constraints.add_line(local_face_dof_indices_fluid[i]);
constraints.add_entry(local_face_dof_indices_fluid[i], column_index, 1.0 / time_step);
constraints.set_inhomogeneity(local_face_dof_indices_fluid[i], - solution_previous_time_step(column_index) / time_step);
}
}

Basically I don't know how to extract the index (or indices) that point at the solid displacement value on the other side of the cell face. Is there a way to do that?

Thank you

Wolfgang Bangerth

unread,
Dec 18, 2024, 9:00:12 PM12/18/24
to dea...@googlegroups.com
On 12/18/24 16:00, Sclah wrote:
>
> for (unsigned int i = 0; i < local_face_dof_indices_fluid.size(); ++i)
> {
> if (fe_fluid->face_system_to_component_index(i).first < dim) // select just
> fluid velocity (the solution is stored as (fluid_velocity + fluid_pressure +
> solid_displacement) with (dim + 1 + dim) components)
> {
> auto column_index = /*todo*/ ; // need to extract the right index to select
> the displacement value at the interface
> constraints.add_line(local_face_dof_indices_fluid[i]);
> constraints.add_entry(local_face_dof_indices_fluid[i], column_index, 1.0 /
> time_step);
> constraints.set_inhomogeneity(local_face_dof_indices_fluid[i], -
> solution_previous_time_step(column_index) / time_step);
> }
> }

Sclah:
if I understand you right, then what you want is to find the displacement DoF
on one side that corresponds to the velocity DoF on the other side? Are you
using the same finite element for the two variables? If so, just as you
extract the vector component with fe_fluid->face_system_to_component_index(),
you can go the other way around with fe_solid->face_component_to_system_index().

Best
WB

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


Sclah

unread,
Dec 19, 2024, 1:52:34 PM12/19/24
to deal.II User Group
Thank you for your help.Yes, I am using the same finite element for the two variables.
I just checked and face_component_to_system_index() doesn't seem to exist inside the library, will component_to_system_index() work just fine?

Wolfgang Bangerth

unread,
Jan 3, 2025, 11:44:54 AM1/3/25
to dea...@googlegroups.com
On 12/19/24 11:52, Sclah wrote:
> Thank you for your help.Yes, I am using the same finite element for the two
> variables.
> I just checked and face_component_to_system_index() doesn't seem to exist
> inside the library, will component_to_system_index() work just fine?

Ah, that's a bummer. When you call
cell->face(f)->get_dof_indices(face_dof_indices)
you get a list of global DoF indices, say
[12, 13, 14, 27, 28, 29]
and you want to find out which of these corresponds to, say, an x-velocity
DoF. fe.system_to_component_index() can tell you that if the list was a list
of DoF indices for the *cell*, but you have a list of DoF indices on the *face*.

Since there is no face_system_to_component_index(), you have two options:
* You write a function FiniteElement::face_system_to_component_index(), using
FiniteElement::system_to_component_index() as an example.

* You find a way to translate from face to cell DoF indices, for example to
know that the 14 above is the face DoF index for shape function 2 (starting
counting at zero) on the cell and that that is, say, shape function 8 on the
cell. Then you can use FiniteElement::system_to_component_index(8) to query
the properties of the 8th shape function. You can make that translation by calling
cell->get_dof_indices(cell_dof_indices)
and then finding where each of the elements of face_dof_indices is located in
cell_dof_indices.

That said, looking this up just now, there *is* a function called
FiniteElement::face_system_to_component_index(). Why can't you use that? There
is not one for the other direction, i.e., there is no
FiniteElement::face_component_to_system_index() but you only need things in
one direction and can get the other direction as the inverse of that mapping.

Best
W.

Sclah

unread,
Jan 10, 2025, 10:16:40 AM1/10/25
to deal.II User Group
Got it, thank you for your help! I will go with the last suggestion using face_system_to_component_index() function and taking the inverse mapping. I didn't think about that, forgive me, I am quite new to dealii library.
Thank you again!
Reply all
Reply to author
Forward
0 new messages