Weak form involving term that depends on unknown solution evaluated at a point

I want to find solution to a Poisson problem, for domain \Omega = (0,1)^2 \subset \mathbb{R}^2:

- \nabla ^2 \mathbf{u} = f \ \ \textrm{in} \ \Omega ,
\mathbf{u} = 0 \ \ \textrm{on} \ \ (0,y),
u_x = c \ \ \textrm{on} \ \ (1,y),

where \mathbf{u} = (u_x,u_y), f=|u(0.5,0.5)| and c>0, c \in\mathbb{R}

Following is the MWE for defining the variational problem:

mesh = UnitSquareMesh(10, 10)
V = FunctionSpace(mesh, "Lagrange", 1)

c = Constant ( 0.1 )
right_fc = CompiledSubDomain("near(x[0],1.0) && on_boundary")               # right face
left_fc = CompiledSubDomain("near(x[0],0.0) && on_boundary")                 # left face
bc_left = DirichletBC(V, Constant((0.0 , 0.0 )), left_fc)
bc_right =  DirichletBC(V.sub(0), c , right_fc)
bc = [ bc_left, bc_right]

# Define variational problem
u = TrialFunction(V)
v = TestFunction(V)
up=u(0.5,0.5)
f = math.sqrt (( up[0] ) **2. +  ( up[1] ) **2.)
a = inner(grad(u), grad(v))*dx
L = f*v*dx 

# Compute solution
u = Function(V)
solve(a == L, u, bc)

I am unsure of the right way of evaluating up=u(0.5,0.5). Could someone please help me with this?

Hi there,

not really your desired solution but for an quick solution I would try to implement a fix point iteration:

  1. Set a start value for f.
  2. Solve the problem
  3. Evaluate f and update its value on the weak form

Do this on a loop until the the value of f converges within a desired tolerance.

Thank you for the reply.

I was rather hoping for find a way of introducing up = u( 0.5, 0.5) in the variational formulation.

I came across this post and was wondering if using Dirac-Delta would work for this problem.

I guess adding u*v*delta*dx, (Where delta is the delta function) to the Left hand side of your problem would be an approximation to having a «source» at the given point.

I was rather thinking of something like (u * delta(0.5 , 0.5) * dx) * v *dx.

then Delta(0.5,0.5) will be 1, and this would be the same as having inner(u,v)*dx

Wouldn’t u * delta(0.5 , 0.5) * dx be the same as \int_{\mathbf{x}_0- \epsilon }^{\mathbf{x}_0+\epsilon} \mathbf{u} (\mathbf{x}) \ \delta (\mathbf{x-\mathbf{x}}_0)d\mathbf{x} = \mathbf{u} (\mathbf{x}_0 ) where \mathbf{x}_0=(0.5 , 0.5) ?

In that case isn’t inner((u * delta(0.5 , 0.5) * dx) , v )*dx the same as inner( u(0.5 , 0.5), v )*dx? Am I missing something?

I probably wasnt clear enough when i wrote the first post if delta(0.5,0.5) is the spatially varying delta function around 0.5,0.5 that is what you would like yes.

1 Like