Hi there, I am starting with the test case on 2D incompressible flow around the cylinder. I went through the tutorial many times but did not figure out where the initial condtions are set.
To be more specific, I would like to do a real simulation of an uniform flow entering the rectangular domain by the inlet boundary. I have already modified the inlet boundary conditions. As for the initial conditions, (v,p) should be zero in the whole domain, and (u) is preferred to be initialized as u_in on the inlet boundary and zero elsewhere. Is that implemented in the tutorial? Otherwise is it possible to set up such an initial condition for (u) ?
As I am not sure what version of FEniCS you are using (legacy FEniCS or DOLFINx), and how you have implemented the flow simulation (direct coupled or via a splitting scheme) I will not provide specific code.
You can use bc.apply(u.vector()) in the case of legacy fenics, and dolfinx.fem.petsc.set_bc(u.vector, [bc]) in the case of DOLFINx to apply a boundary condition to a function.
I am using DOLFINx. petsc.set_bc works well. The simulation seem to diverge rapidly, but it’s an irrelevant topic to this one…Thank you anyway.
When I tried to transform the velocity from function to numpy arrays, I get strange array shapes. Indeed, u_.x.array.shape returns (25478,), but p_.x.array.shape and mesh.geometry.x.shape return (3252, ) and (3252, 3). I recognize 3252 as the number of nodes in my mesh. The returned value about velocity seem to be caused by the CG2 function space for velocity. Do you know how can I get the numpy array of velocity fields ?
What do you mean by the numpy array of the velocity field?
When you use a CG-2 function space, you have 6 degrees of freedom (see: DefElement)
Thus it will not match the number of vertices of your mesh.
Also note that each component (x,y,z) of the velocity is a separate degree of freedom.
The array u_.x.array matches the coordinates of (x=V_space.tabulate_dof_coordinates()), where the i-th row of x corresponds to the u_.x.array[2*i:2*(i+1)] entries of u_ (in the 2D setting).
I try to obtain the velocity values on mesh nodes, like what I get by calling p_.x.array for the pressure. V_space.tabulate_dof_coordinates() looks like a refined version of the original mesh.geometry.x, but I do not understand the transformation in between. The correspondence describe in DefElement doesn’t match the reality, at least the first 10 rows of V_space.tabulate_dof_coordinates() and mesh.geometry.x.
Doesn’t there exist a clean way to obtain velocity values on the mesh nodes?
The point is that there is no one-to-one correspondence between the dof coordinates and the mesh vertices, as illustrated by the plots on DefElement. A second order space has dofs on the facets and vertices. If you want to have only values at the nodes, you could interpolate your CG-2 solution into CG-1.
See for instance: Mapping of results to mesh nodes, debugging of dolfinx implentation of shell formulation - #6 by dokken and the previous posts in that thread