Acoustical problems with symmetry in ngsolve

More
4 years 6 months ago #2700 by gcdiwan
Dear NGSolve developers;

I am trying to understand how NGSolve solves the problems with boundary conditions defined on surface that is symmetric about a coordinate plane. When solving
[tex]-\Delta u -k^2 u = 0 \quad {\rm{\in}} \quad \Omega \\
\frac{\partial u}{\partial n} = 1.0 \quad {\rm{on}} \quad \Gamma_1\\
\frac{\partial u}{\partial n} -\beta u = 0 \quad {\rm{on}} \quad \Gamma_2[/tex]
with \Gamma_1 (symmetric_srf) and \Gamma_2 (outer) nonintersecting, \Omega a sphere of given dimensions, and \Gamma_1 is symmetric about one of the coordinate planes. Taking advantage of the symmetry of \Gamma_1 when the problem is divided in half, I see the solution is scaled by a factor of 2 compared to the original problem. This is a bit surprising because in the original problem the pressure field is symmetric, i.e., du/dn = 0 on the symmetry plane. Hence I should have been able to take any subset of the full problem, and apply either Dirichlet or Neumann boundary conditions, taken from the full problem, at the boundary of the subset, and be able to compute the solution correctly, using the boundary conditions and sources internal to the subset. Note that when the full model is sliced into half to take advantage of the symmetry, i don't specify any condition on the newly created surface (thus a natural bc is imposed there which i believe is already from the subset of the original problem). Could you comment on what i could be doing wrong here?
Below the relevant snippet
Code:
mesh = Mesh(ngmesh) fes = H1(mesh, complex=True, order=2) u = fes.TrialFunction() v = fes.TestFunction() g = 1; beta = 1j*waveno a = BilinearForm (fes, symmetric=True) a += SymbolicBFI (grad(u)*grad(v) ) a += SymbolicBFI (-waveno*waveno*u*v) a += SymbolicBFI(beta*u*v, definedon=mesh.Boundaries("outer")) f = LinearForm (fes) f += v * g * ds("symmetric_srf") gfu = GridFunction (fes) a.Assemble() f.Assemble()
More
4 years 6 months ago #2701 by christopher
At first glance I do not see any problems with that. Could you post a minimal example?
Best
More
4 years 6 months ago #2709 by gcdiwan
Hi Christopher,

MWE:
I tried to recreate the pressure doubling issue with a different geometry -
half and full model for sphere (R1=1) in sphere (R2=2) where the interior sphere surface is excited with neumann bc. I tried both an impedance and PML condition. I don't see pressure doubling in the MWE which is reassuring.

Real Geometry where pressure doubling happens:
Unfortunately, I can't post the original geometry STEP file here but the
GMSH mesh for the actual geometry is attached - realGeomHalfmodel.txt.
The attached code.py is practically what i used for the real geometry problem (except the bit where python calls gmsh and meshes using .geo file) where i get the pressure doubling. So i am guessing i must be doing something wrong when i create the mesh via the .geo file. Plane 6 is symmetrical about the vertical plane passing through it and is excited with neumann bc. In the full model, plane 6 would have been a full circle and the model would have contained the elements from the other half.
More
4 years 6 months ago - 4 years 6 months ago #2710 by christopher
If you can send me the step file at clackner@cerbsim.com ? The attached mesh seems to have no surface elements.
Is there a specific reason you use gmsh for the meshing? If you use netgen you can have higher order curved boundaries which would be better for absorbing bc on a sphere...
My first guess is that something in the mesh is not correct that may cause the factor.
Best
Last edit: 4 years 6 months ago by christopher.
More
4 years 6 months ago #2711 by joachim
Your observation is correct:

If you apply the source at the symmetry plane of the full model, the weak form leads to the interface condition

[du/dn] = g

i.e.
du/dn_left + du/dn_right = g,

i.e. by symmetry
du/dn_left = du/dn_right = g/2

while for the half problem the source term is
du/dn = g

in physical words:
if you apply the source term for the full problem, half of the energy is going to the left, and half to the right.
if you do the half problem, the full energy is going to this side.

Joachim
More
4 years 6 months ago #2712 by gcdiwan
The reason for using gmsh is ease with which geo files can be edited and also ability to assign various physical tags to the surfaces/lines/volumes etc.
Time to create page: 0.111 seconds