1D problem with Neumann Boundary Conditions


98
views
0
3 months ago by

I am working on applying Neumann Boundary Conditions to a 1D problem with no luck. A minimal working example is below. Any tips on what I am doing wrong would be greatly appreciated.

from fenics import *
import numpy as np

n = 200

mesh = UnitIntervalMesh(n)
V = FunctionSpace(mesh, 'CG', 1)
v = TestFunction(V)
rho = Function(V)

class LeftEdge(SubDomain):
	def inside(self, x, on_boundary):
		return on_boundary and near(x[0], 0)

class RightEdge(SubDomain):
	def inside(self, x, on_boundary):
		return on_boundary and near(x[0], 1)

boundaryLeft = LeftEdge()
boundaryRight = RightEdge()

boundaries = FacetFunction("size_t", mesh)
boundaries.set_all(0)
boundaryLeft.mark(boundaries,1)
boundaryRight.mark(boundaries,2)

ds = ds(subdomain_data = boundaries)

G = Constant(1)

F = rho * v * dx \
  + v * G * ds(1) \
  - v * G * ds(2)

solve(F == 0, rho)

rho = rho.vector().get_local()

print('Resulting BC', np.diff(rho, 1./n)[[0,-1]], 'Applied BC', G.values()[0])
Community: FEniCS Project

1 Answer


3
3 months ago by
Your variational form F doesn't involve any derivatives.  You cannot impose boundary conditions without a PDE.
That's right, from you weak form, it seems you want to solve \(u = 0\) with Dirichlet boundary conditions. That's why you get a highly oscillatory solution – this is more obvious if you use a discontinuous function space. Also, np.diff() gives you the difference, not the derivative.
written 3 months ago by Adam Janecka  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »