I've found what was wrong in this piece of code. It turned out that 
index.first < dim does not guarantee this dof is a velocity dof. Changing 
it to* if (i < fe.base_element(0).dofs_per_cell && index.first < dim)* 
solves the problem.

But there are two things that I can not understand:
1. using i < fe.base_element(0).dofs_per_cell to make sure the current dof 
is a velocity dof which I want to constrain, but this is obviously not good 
practice as it assumes the velocity dofs are always in front of pressure 
dofs. What is the correct way?
2.  index.first should return the component number: if the dof is a 
velocity dof then it would be in range [0, dim), if it is a pressure dof 
then it should be 0. However, based on my test index.first can be greater 
than dim, why? Then how to get the velocity component number (0 for x, 1 
for y)?

*    if (index.first < dim) // Is this the correct method to tell if the 
> dof i is a velocity dof rather than a pressure dof??*
>     { 
>       Point<dim> mapped_point = mapping.transform_unit_to_real_cell(cell, 
> unit_points[i]); // get the real coordinates
>       Vector<double> value = f(mapped_point); // Calculate (interpolate) 
> the constrained value which is a *vector*
>     
>       // set the inhomogeneity constraints on the velocity
>       *auto line = dof_indices[i]; // I think this is the global dof 
> number corresponding to dof i*
>       constraints.add_line(line);
>       *constraints.set_inhomogeneity(**line, value[index.first]); // the 
> component number should be the correct velocity component to constrain *
>     }
>

Thanks
Jie 

-- 
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+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Reply via email to