Forum Message

 

 

We have moved the forum to https://forum.ngsolve.org . This is an archived version of the topics until 05/05/23. All the topics were moved to the new forum and conversations can be continued there. This forum is just kept as legacy to not invalidate old links. If you want to continue a conversation just look for the topic in the new forum.

Notice

The forum is in read only mode.

Boundary definitions and Dirichlet BCs on specific DOFs for vector valued funcs

More
4 years 8 months ago #1813 by anandes
Hello,

I'm new to Netgen/NGSolve. I'm trying to solve a few simple elasticity problems to start with. I have a few questions in this regard.

1. How do we define boundaries which are only a part of an edge in 2D or face in 3D (for example: not the entire bottom, top, etc. portion of the boundary)?
1a Can this be done in the Python interface based on the geometric location of the portion of the boundary surfaces where the Dirichlet BCs need to be imposed?
2. Is it possible to impose Dirichlet BCs at a specific node? If so, how is this done?
3. How do we impose Dirichlet BCs on a specific DOF of a given node?

Thank you,
Anand
More
4 years 8 months ago #1816 by christopher
Hi Anand,
1. you can split the boundary in the geometry into multiple parts, each using one segment and then defining the bc on the subpart of the edge.
2&3. Yes you can query the dirichlet dofs of a finite element space and set them to not free if they should be dirichlet dofs
Code:
freedofs = fes.FreeDofs() freedofs.Clear(dofnr)
You can iterate over the mesh topology, as a reference see here:

ngsolve.org/docu/latest/i-tutorials/unit...gy/meshtopology.html

and for example set all edge and vertex dofs of edges to dirichlet where both vertices have an x-coordinate smaller than 0:
Code:
for edge in mesh.edges: if all(mesh[v].point[0] < 0 for v in edge.vertices): for vertex in edge.vertices: for dof in fes.GetDofNrs(vertex): freedofs.Clear(dof) for dof in fes.GetDofNrs(edge): freedofs.Clear(dof)
Note, that mesh.edges only gives you NodeIDs, you have to get the node by inserting them into the mesh and then you can get the point (and then the x value of that) from a vertex node.

Best
Christopher
More
4 years 8 months ago #1822 by anandes
Hi Christopher,

Thanks for your detailed clarifications. I'll try this and let you know if I face issues.

Thank you,
Anand
More
4 years 8 months ago #1826 by anandes
Hello Christopher,

I'm following up on my earlier questions. I've attached a file with a 2D mesh (SEC_2D_lin_tri.vol) and the python script (try2.py) to set up a simple static elasticity problem in NGSolve. The solution is a 2D vector field with components, say, ux & uy. The Dirichlet BCs to be applied are:
a) uy=0 on Ysym (x>=0 and y=0)
b) ux=0 at the bottom right corner (x=10, y=0)

The top surface (Pres_load) has a non-zero Neumann condition. This is a well-posed problem and is a standard test problem in linear elasticity.

When I use the function space for the trial and test functions and set the flag dirichlet="Ysym", both degrees-of-freedom (ux and uy) are enforced on Ysym - this is not what I want. So I do not set the dirichlet flag while declaring the funtion space. Instead, I use lines 34 to 55 in the python script to define the Dirichlet BCs similar to your suggestion. When I print the FreeDofs() bit array, all the bits corresponding to the Dirichlet BCs are '0', and I verified that the correct DOFs are selected. However, the system matrix seems to be singular and the solution is NaN everywhere.

Can you let me know what I'm getting wrong here? Is there a cleaner way to implement BCs on individual DOFs (something similar to line 92)?

Another unrelated issue: When I read this mesh and try a uniform refinement, mesh.Refine(), the code throws a segmentation fault.

Thanks,
Anand

File Attachment:

File Name: SEC_2D_lin_tri.vol
File Size:14 KB

File Attachment:

File Name: try2.py
File Size:3 KB
More
4 years 8 months ago - 4 years 8 months ago #1831 by christopher
The VectorH1 space already has as many components as the mesh dimension. If you set dim there you will get dim * mesh.dim dimensions, this is not what you want.
Code:
fes = VectorH1(mesh, order=1)

You can set your boundary conditions this way:
Code:
fesx, fesy = fes.components for el in fesx.Elements(BND): if el.mat == "Ysym": for dof in el.dofs: freedofs.Clear(fesx.ndof + dof) # offset of x-dofs # find bottom right vertex for v in mesh.vertices: if abs(v.point[0] - 0.5*SB) + abs(v.point[1]) < tol: for dof in fesx.GetDofNrs(v): freedofs.Clear(dof)

Best
Christopher
Last edit: 4 years 8 months ago by christopher.
More
4 years 8 months ago #1835 by anandes
Hello Christopher,

Thanks for your response. I'm now able to apply Dirichlet BCs on specific DOFs at a given node.

If I try to use your method with the 'dim' parameter instead of VectorH1 - NGSolve throws this error " 'components' is available only for product spaces". What is the difference between setting VectorH1 and setting the 'dim' parameter? If it is too much to explain on this forum can you please point me to any NGSolve documentation if available?

Regarding mesh refinement for the previously attached .vol mesh file: as I mentioned earlier, NGSolve throws a segmentation fault when I try mesh.Refine(). Can you please let me know what the issue could be?

Thanks,
Anand
Time to create page: 0.179 seconds