Implementation of Neumann boundary condition in 3d

[HIMAL MAGAR]

Greetings,

I am trying to apply Neumann boundary condition(pressure) on inner wall of a cylinder. So far, I have followed Lecture 21.55 to apply this condition in 2d (in step-6 and step-8) and succeeded. But when it came to 3d case, I am having hard time understanding how to implement pressure which is always normal to the faces of inner wall.

Thank you.

Himal.

[Wolfgang Bangerth]

Himal:
the first step in these cases is always to write down the weak formulation of
the equation. Can you explain what the equation is and where Neumann boundary
conditions enter? There shouldn't really be a difference between 2d and 3d,
you will generally just end up with a term of the form
(p, n.phi)_{Gamma_N}
where p is the prescribed pressure, n the normal vector, and phi a
(vector-valued) test function, all integrated over that part of the boundary
Gamma_N where you impose Neumann boundary conditions. The form of this
integral does not depend on whether you are in 2d or 3d.

Best
W>

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

[HIMAL MAGAR]

Dear Sir,

The equation is Poisson's equation with mixed boundary conditions. The domain is a cylinder shell. I am trying to apply Homogeneous Boundary Condition on the bottom face as done in step-18 and Neumann Boundary Condition on the inner surface of the cylinder. I have already derived a weak formulation of the equation and come up with the term you have mentioned. When applying in the code, do I have to construct a separate template for "pressure"? or can i just put the value explicitly in the cell_rhs component for Neumann boundary? And if I do the latter, what will be the direction of the pressure relative to the face?

Thank you.

Himal

--
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 a topic in the Google Groups "deal.II User Group" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/dealii/_4eeywhhxk4/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dealii+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/eb88ce6b-90bb-8dcc-c65b-6556f9f2b17c%40colostate.edu.

[HIMAL MAGAR]

Dear Sir, 

So far I have used step-17 and step -18 to apply body force and zero Dirichlet b.c.'s.  I have attached the result below:-).

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?

Thanking you.

Himal.

[Wolfgang Bangerth]

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.

[HIMAL MAGAR]

Dear SIr,

Sorry for the late response. Yes, that makes sense. So,
I have been working on this case and faced different problems.

I have used step-20 and your suggestion to implement the neumann b.c. used the code as follows:

The result was output in vector format following your Lecture 21 using step-22 which is as follows:

Here the scale was set to 130.

The same simulation was done in SolidWORKS and following result was obtained:

Scale was the same.

I am extremely happy about the result :-).

For now, I want to output norm_stress in the results. In step-18, it is said that the process is simple, however, i have found that difficult to incorporate in my code. How can I approach this, Sir?

Thank you.

Himal.

--
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 a topic in the Google Groups "deal.II User Group" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/dealii/_4eeywhhxk4/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dealii+un...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/8fa32629-8350-b88c-7d42-6131a4fba57d%40colostate.edu.

[Wolfgang Bangerth]

On 11/21/22 03:04, HIMAL MAGAR wrote:
> 0-02-03-ecf50731d20a4fed58e332ff9a1702a1e061032d7400f34e6ee11a2b55bd6aea_2073f8a6123d5b.png
> Scale was the same.
> I am extremely happy about the result :-).

Very nice!

> For now, I want to output norm_stress in the results. In step-18, it is said
> that the process is simple, however, i have found that difficult to
> incorporate in my code. How can I approach this, Sir?

Ask more specific questions. "I have found that difficult" is not a
description of the problems you have faced, and it's difficult to tell you
what exactly you need to do to make it work.

Separately, you might want to take a look at the DataPostprocessorScalar
class. step-29 provides an example of how one can output the norm of a
solution -- in that case, the norm of the solution, not its gradient.