Is it possible to do Vibroacoustics ?

More precisely i want to know if there is any pitfall for the implementation of something like seen in this elmer fem video.



It would be great if i can use ngsolve for this task

Hi Kai,

I guess you want to solve it in two steps ?
Step 1 eigenvalue problem for the tuning forg.
Step 2 solve the acoustic problem

For step 1 you may watch this video from a student's seminar:


For step 2 you may want to use PML for the radiation condition, see
docu.ngsolve.org/latest/i-tutorials/unit-1.7-helmholtz/pml.html

Best, Joachim

Ok, i have played around with some examples where the wave equation is used. I finally end up with a code that have the boundaries of a rectangle as source for the wave equation.
[/code]

I now want to use my computed eigen modes for the source geometry, that i compute as
[/code]


How can i achieve this ?

I further have some problems for the 3D version of all this, i do not really understand how to asign boundary conditions to specific parts of a herachical mesh. My code so far to generate the acoustic structure mesh is
[/code] So i would be very happy if some one could give me some hints to put this all together

 

hm ... my code is not visible in firefox ...

ok -repost for check if i had made a formating mistake
Ok, i have played around with some examples where the wave equation is used. I finally end up with a code that have the boundaries of a rectangle as source for the wave equation.

from ngsolve import *
from netgen.geom2d import SplineGeometry

# Geometry
geo = SplineGeometry()
geo.AddCircle((0.0, 0.0), 0.8, bc="outer")
geo.AddRectangle((-0.05, -0.025), (0.025, 0.275),
leftdomain=0, rightdomain=1, bc="scat")
mesh = Mesh(geo.GenerateMesh(maxh=0.1))

fes = H1(mesh, order=5, complex=True)
u, v = fes.TnT()

# Wavenumber & source
omega = 15
pulse = 5e4*exp(-(40**2)*((x-0.5)*(x-0.5) + (y-0.5)*(y-0.5)))

# Forms
a = BilinearForm(fes)
a += grad(u)*grad(v)*dx - omega**2*u*v*dx
a += -omega*1j*u*v * ds("outer")
a.Assemble()

f = LinearForm(fes)
f += 1 * v * ds("scat")
f.Assemble();

gfu = GridFunction(fes, name="u")
gfu.vec.data = a.mat.Inverse() * f.vec
Draw(gfu)


I now want to use my computed eigen modes for the source geometry, that i compute as

a, b, fes = build_elasticity_system(mesh, steel, dirichlet)
 u = ngsolve.GridFunction(fes, multidim=num)
lams = ngsolve.ArnoldiSolver(a.mat, b.mat, fes.FreeDofs(), u.vecs, shift)


How can i achieve this ?

I further have some problems for the 3D version of all this, i do not really understand how to asign boundary conditions to specific parts of a herachical mesh. My code so far to generate the acoustic structure mesh is

import ngsolve

from ngsolve import *

from netgen.stl import *
from netgen.csg import *
from netgen.meshing import *

geo2 = STLGeometry("../test_data/tuningfork.stl")
m2 = geo2.GenerateMesh()

centerx = 0
centery = 0
centerz = 0

count = 0

for e in m2.Elements2D():
    for v in e.vertices:
            centerx = centerx + m2[v][0]
            centery = centery + m2[v][1]
            centerz = centerz + m2[v][2]
            count = count + 1

centerx = centerx / count
centery = centery / count
centerz = centerz / count

R = 0
for e in m2.Elements2D():
    for v in e.vertices:
            x = -centerx + m2[v][0]
            y = -centery + m2[v][1]
            z = -centerz + m2[v][2]
            tmp = sqrt(x*x + y*y + z*z)
            if tmp > R:
                R = tmp

print("R : " + str(R))

geo1 = CSGeometry()
# Set the mesh size on the sphere surface to 0.1
sphere = Sphere(Pnt(centerx, centery, centerz), 2.5*R)
geo1.Add(sphere)
m1 = geo1.GenerateMesh()

print(m2)

# create an empty mesh
mesh = netgen.meshing.Mesh()

# a face-descriptor stores properties associated with a set of surface elements
# bc .. boundary condition marker,
# domin/domout .. domain-number in front/back of surface elements (0 = void),
# surfnr .. number of the surface described by the face-descriptor

fd_outside = mesh.Add(FaceDescriptor(bc=1,domin=1,surfnr=1))
fd_inside = mesh.Add(FaceDescriptor(bc=2,domin=2,domout=1,surfnr=2))
# copy all boundary points from first mesh to new mesh.
# pmap1 maps point-numbers from old to new mesh

pmap1 = {}
for e in m1.Elements2D():
    for v in e.vertices:
        if v not in pmap1:
            pmap1[v] = mesh.Add (m1[v])


# copy surface elements from first mesh to new mesh
# we have to map point-numbers:

for e in m1.Elements2D():
    mesh.Add(Element2D(fd_outside, [pmap1[v] for v in e.vertices]))

# same for the second mesh:
pmap2 = {}
for e in m2.Elements2D():
    for v in e.vertices:
        if v not in pmap2:
            pmap2[v] = mesh.Add(m2[v])

n = len(pmap2)

for e in m2.Elements2D():
    mesh.Add(Element2D(fd_inside, [pmap2[v] for v in e.vertices]))


mesh.GenerateVolumeMesh(maxh=10.0)
mesh.Save("acousticmesh.vol")

mesh = ngsolve.Mesh(mesh)

So i would be very happy if some one could give me some hints to put this all together



 

looks like the [code] tag does not work