| 2013-07-30 15:20:57 |
Peter Kekenes-Huskey |
description |
Hi,
I modified the
# Heavily borrowed from
# http://fenicsproject.org/pub/book/tutorial/stationary/poisson/mat2_p2D.py
"""
FEniCS tutorial demo: Poisson problem in 2D with 2 materials.
"""
Hi,
I tried to apply the protocol for marking cells using the based on the poisson problem with two materials example (http://fenicsproject.org/pub/book/tutorial/stationary/poisson/mat2_p2D.py) for a 2D mesh I created with gmsh (attached). It appears that the classes for marking the subdomains do not work for my mesh (several elements are missed); unfortunately it is not obvious to me why this is failing. Is there something incorrect about the mesh I am using (it runs fine with other dolfin PDE solves) or the code I have below?
Thanks for your help
pete
from dolfin import *
import sys, math, numpy
import numpy as np
class Omega0(SubDomain):
def inside(self, x, on_boundary):
return True if x[1] <= 0.5 else False
class Omega1(SubDomain):
def inside(self, x, on_boundary):
return True if x[1] >= 0.5 else False
if(1):
nx = 4; ny = 6
mesh = UnitSquare(nx, ny)
#mesh = Mesh('nohole.xml')
mesh = Mesh('hole.xml')
# Define a MeshFunction over two subdomains
subdomains = MeshFunction('uint', mesh, 2)
# Mark subdomains with numbers 0 and 1
subdomain0 = Omega0()
subdomain0.mark(subdomains, 0)
subdomain1 = Omega1()
subdomain1.mark(subdomains, 1)
# Discontinuous basis is used for boundary condition
# LAter, continuous basis used for PDE solution
V0 = FunctionSpace(mesh, 'DG', 0)
k = Function(V0)
print 'mesh:', mesh
print 'subdomains:', subdomains
print 'k:', k
import numpy
help = numpy.asarray(subdomains.array(), dtype=numpy.int32)
#print help
print "bad ", np.where(help>1)
k.vector()[:] = numpy.choose(help, k_values) |
Hi,
I tried to apply the protocol for marking cells using the based on the poisson problem with two materials example (http://fenicsproject.org/pub/book/tutorial/stationary/poisson/mat2_p2D.py) for a 2D mesh I created with gmsh (attached). It appears that the classes for marking the subdomains do not work for my mesh (several elements are missed); unfortunately it is not obvious to me why this is failing. Is there something incorrect about the mesh I am using (it runs fine with other dolfin PDE solves) or the code I have below?
Thanks for your help
pete
from dolfin import *
import sys, math, numpy
import numpy as np
class Omega0(SubDomain):
def inside(self, x, on_boundary):
return True if x[1] <= 0.5 else False
class Omega1(SubDomain):
def inside(self, x, on_boundary):
return True if x[1] >= 0.5 else False
if(1):
nx = 4; ny = 6
mesh = UnitSquare(nx, ny)
#mesh = Mesh('nohole.xml')
mesh = Mesh('hole.xml')
# Define a MeshFunction over two subdomains
subdomains = MeshFunction('uint', mesh, 2)
# Mark subdomains with numbers 0 and 1
subdomain0 = Omega0()
subdomain0.mark(subdomains, 0)
subdomain1 = Omega1()
subdomain1.mark(subdomains, 1)
# Discontinuous basis is used for boundary condition
# LAter, continuous basis used for PDE solution
V0 = FunctionSpace(mesh, 'DG', 0)
k = Function(V0)
print 'mesh:', mesh
print 'subdomains:', subdomains
print 'k:', k
import numpy
help = numpy.asarray(subdomains.array(), dtype=numpy.int32)
#print help
print "bad ", np.where(help>1)
k.vector()[:] = numpy.choose(help, k_values) |
|
