On 11/17/22 09:36, HIMAL MAGAR wrote:
>
> I am now trying to implement Neumann b.c.'s on inner wall.
> So this is what i am thinking:
> #Create a template such that it returns a vector with constant magnitude
> and direction indicated by unit vector along radial direction. Does it
> make sense?
This is difficult because the function you have to return depends not
only on the pressure but also on the normal vector -- but the latter
depends on where you are on the geometry.
In the end, when you compute the integral for the boundary term in the
case of elasticity, you will have to compute something like
(\vec g, \vec \phi)_{Gamma_N}
which you have to do via quadrature. In your case, you have
\vec g = p \vec n
So you would do something like this:
for (const auto &cell : dof_handler.active_cell_iterators())
for (const auto &face : cell->face_iterators())
if (face->at_boundary() && is inner face of your cylinder)
{
fe_face_values.reinit(cell
for (q=0....n_q_points)
{
const Tensor<1,dim> n = fe_face_values.normal_vector(q);
const double p = ...whatever the pressure is here...;
const double phi_i
= fe_face_values[displacement].shape_value(i,q);
rhs_values[i] += (n*p*phi_i) * fe_face_values.JxW(q);
}
... copy local to global ...
Does this make sense?
Best
W.