Dear Martin,
> This procedure does not give me a consistens right hand side as the
> part g in rhs=[b;g]'; does not satisfy
> (g,1)=0;
I've just had a look at the pressure right hand side from the step-22
and it is, as in your case, not equal to zero. And it is not a problem
with component_mask or the way you define your BoundaryValues function,
at least as far as I can judge.
The reason for having some non-zero entries in the vector g, which is
the right hand side of
div(u) = 0,
in case that we don't have any boundary value on the pressure, is
because of what you're doing with the velocity:
Some velocity (!) degrees of freedom are fixed to a non-zero value and
they cause some entries. To understand this, let's have a look at a
(weird) 1D equation of this form with conditions u(0) = u(1) = 1. The
discrete equation for linear FE's is basically:
0.5 (u_1 - u_0) = 0,
0.5 (u_{i+1} - u_{i-1}) = 0 at interior nodes, i=1,...,N-1,
0.5 (u_N - u_{N-1}) = 0.
The boundary conditions set u_0 = 1, u_N = 1.
The ConstraintMatrix then substitutes u_0 and u_N from the above
difference equation and puts the known value to the right hand side, so
the equation the you get after the condense command actually reads:
0.5 u_1 = 0.5,
0.5 u_2 = 0.5,
0.5 (u_{i+1} - u_{i-1}) = 0, i=2,...,N-2,
0.5 u_{N-2} = 0.5,
0.5 u_{N-1} = 0.5.
So the right hand side is non-zero!
Exactly the same happens in 2D/3D when there are some velocity
components that are fixed to a non-zero boundary value. So I think that
the non-zero vector g that you encounter doesn't look like an error,
since there is some natural explanation...
Best,
Martin
_______________________________________________
dealii mailing list http://poisson.dealii.org/mailman/listinfo/dealii
