Assembling of RHS for nonlinear Helmholtz eq. (Thermoacoustics)

Hello, I’m implementing a solver for the eigenvalue problem of the nonlinear Helmholtz eq. for thermoacoustics:
\nabla (c^2 \nabla p)+\omega^2 p = FTF(\omega)\nabla p|_{ref} \cdot n_{ref}
where FTF is some Flame function (in my case n-tau model) and n_{ref} is some unit vector.
I’ve had some problems when trying to assemble the RHS. General speaking, the matrix should look like:
F=\int_{flame} \phi_j dV \nabla \phi_i|_{x_{ref}}\cdot n_{ref} =qg^T

in other words:
F=\begin{cases} 0 & x \not\in \Omega_{ref} \\ \int_{flame} \phi_j dV\nabla \phi_i|_{x_{ref}}\cdot n_{ref} & x \in \Omega_{ref} \\ \end{cases}
or q:= as the standard load vector & g defined as :
g=\begin{cases} 0 & x \not\in \Omega_{ref} \\ \nabla \phi_i|_{x_{ref}}\cdot n_{ref} & x \in \Omega_{ref} \\ \end{cases}
with \Omega_{ref}\subset\Omega.
As a result it is a very sparse matrix and it doesn’t feature a diagonal structure, but i has a few entries along some columns.
I’ve tried to define a new set of test functions over this ref. domain only (as a subdomain) or to implement some type of Dirac delta function, but I haven’t had success so far.
maybe you have some ideas how to solve this problem.
Thanks in advance.
Dolfin v. 2019.1.0

I would strongly suggest adding the code that you have used so far, and explain what error message you are obtaining, as this will increase the likelihood of anyone being able to help you.

Thank you for the advice. Here is some example code, where i try to define a second set of test functions with support only in \Omega_{ref} . I do not obtain any error but the matrix F is not what I’m expecting.

As a result it is a very sparse matrix and it doesn’t feature a diagonal structure, but i has a few entries along some columns.

mesh = UnitSquareMesh(4,4) #Mesh definition

class Sub_domain(SubDomain):     #Subdomain definition as [0.4-eps, 0.4+eps]^2
    def inside(self, x, on_boundary):
        eps=1e3
        return between(x[0], (0.4+eps, 0.4-eps)) and between(x[1],(0.4+eps, 0.4-eps))
subdomain=Sub_domain

V = FunctionSpace(mesh, "CG", 1)
W = FunctionSpace(mesh, "CG", 1,subdomain)

u = TrialFunction(V)             #Trial function for general & local doamin
v = TestFunction(V)              #Test function for the general Domain
w = TestFunction(W)              #Test function for the local Domain
  
f = inner(grad(u)[0],w)*dx       #RHS
F = PETScMatrix()   
assemble(f, tensor=F)

a = inner(grad(u),grad(v))*dx   #Stiffness Matrix
A = PETScMatrix()
assemble(a, tensor=A)

Hello,
I have narrowed down my question:
Is It possible to assemble a vector with trial functions \phi_i without the integration. I would like to do this in order to assemble the matrix F via F=(\int_{flame} \phi_j dV )(\nabla \phi_i|_{x_{ref}}\cdot n_{ref}) =qg^T. This would be helpful for the vector g defined as:

I have tried with:

from fenics import *
mesh = UnitSquareMesh(4,4) #Mesh definition
V = FunctionSpace(mesh, "CG", 1)
u = TrialFunction(V)             #Trial function for general & local doamin
f = grad(u)[0]                   # g vector   
F = assemble(f)

Error:


  File "/Library/code/test.py", line 191, in <module>
    F = assemble(f)

  File "/opt/anaconda3/envs/fenicsproject2019_spyder/lib/python3.9/site-packages/dolfin/fem/assembling.py", line 198, in assemble
    dolfin_form = _create_dolfin_form(form, form_compiler_parameters)

  File "/opt/anaconda3/envs/fenicsproject2019_spyder/lib/python3.9/site-packages/dolfin/fem/assembling.py", line 60, in _create_dolfin_form
    raise TypeError("Invalid form type %s" % (type(form),))

TypeError: Invalid form type <class 'ufl.indexed.Indexed'>`Preformatted text`

Assemble is an integration command.
If you want to evaluate a trial function at a given point, see for instance:
https://bitbucket.org/fenics-project/dolfin/src/43642bad27866a5bf4e8a117c87c0f6ba777b196/python/test/unit/fem/test_manifolds.py?at=master#test_manifolds.py-235,243
or
Evaluate_basis() gives strange results - #2 by kamensky (and the original post in the linked thread)

I was try to assemble the matrix the best way possible. At the end I used basic matrix operations to modify the assembled matrix.

Hi @reckx_sanchez,

Can you please share a minimal code for matrix F?