Hi! I am trying to create radially symmetric sources in R^2 for the purpose of building a (more) well-posed inverse problem. My first idea was to create a 1D fes, simulating functions that depend only on the radius r, then extend these to a 2D fes, which I’d use to solve PDEs and so on.
My first simple idea was
from netgen.geom2d import SplineGeometry
from ngsolve.meshes import Make1DMesh
q = GridFunction(fes1)
q.Set(1-x**2)
f = GridFunction(fes2)
f.Set(q)
Well, this obviously fails in that f becomes constant in y, but it has the correct behaviour in x. Is there a way to extend q to a radially symmetric f, or a better way to create a “space” of radially symmetric sources? To be clear, just specifying the sources as radially symmetric CFs will not do, since I will be solving inverse problems into the space of radially symmetric sources.
This does not seem consistent with the function 1-x^2-y^2. However, if I set order = 2 in the 1d fes, I get the below, which does appear to be very close to the exact function 1-x^2-y^2. Is there a clear reason for this?
As a side note, if this does work with order = 2, then I have the mapping from L^2 over 1d to radially symmetric sources in L^2 over 2d. The adjoint operation (which I will need) from L^2 over 2d to L^2 over 1d would be something like a polar coordinate integral over angle only. How would one go about to implement that?
Sorry my bad, sent you code that I didn’t really test. yes with lower order it looks fine. I’m not yet quite sure why high order is so instable here, will investigate. As mentioned, this space is relatively new and not widely adopted/tested, so handle with care
globalinterface space creates a global basis on parameterspace [0,1] using the given mapping to parameter space. You can visualize the basis functions of the space
Draw(q, mesh, "q")
for i in range(len(q.vec)):
q.vec[:] = 0
q.vec[i] = 1
Redraw()
input("next")
you can assemble mass matrix and rhs on the global space and use this to interpolate into the space by solving that system.
works with the symbolic forms:
fes1 = comp.GlobalInterfaceSpace(mesh, mapping=1 - (x**2 + y**2)/r**2, order=3)
u,v = fes1.TnT()
m = BilinearForm(u * v * dx)
m.Assemble()
f = LinearForm(func_on_circle * v * dx)
f.Assemble()
gf = GridFunction(fes1)
gf.vec.data = m.mat.Inverse() * f.vec
gf would then be the average of func_on_circle over angle.