4
$\begingroup$

I am confused most of the time with the application of boundary conditions in Fenics. Suppose I have an xml mesh for a rectangular beam and want to apply Dirichlet BC ($u =0$) on the left side of the beam in my mesh file. Can anyone tell me how can I achieve this?

I know another way of doing it by defining subdomains like here:

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

But for complicated geometries, I am sure this method is not useful.

Thanks.

$\endgroup$
8
$\begingroup$

That depends how did you generated the file. For example your mesh generator can produce some cell or facet markers according to its input. Then dolfin-convert script may or may not succeed converting all these markers to XML. They're then stored within XML mesh file or in separate XML. You can then use them as FacetFunction

# if stored within mesh
facet_domains = mesh.domains().facet_domains(mesh) # older versions of DOLFIN
facet_domains = mesh.domains().facet_domains() # newer versions of DOLFIN

# if stored separately
facet_domains = MeshFunction('uint', mesh, file) # older versions of DOLFIN
facet_domains = MeshFunction('size_t', mesh, file) # newer versions of DOLFIN

# now define homogenous Dirichlet BC where facet_domains is equal to 42
bc42 = DirichletBC(V, 0.0, facet_domains, 42)

# now integrate over interior facets where facet_domains equals 43
dS = Measure('dS')[facet_domains]
dS43 = dS(43)
integral43 = foo*dS43

If your mesh generator does not supply you with markers you need, you can produce first these markers using EntityIterators, for example

facet_domains = FacetFunction('size_t', mesh, 0)
for c in cells(mesh):
    # here some conditions on c
    for f in facets(c):
        # here some conditions on f
        facet_domains[f] = some_computed_value

Then you can use these facet_domains as suggested above.


Application of DirichletBC at one particular vertex

Defining DirichletBC in one vertex is not possible through MeshFunctions as DirichletBC only supports FacetFunctions. But you can take advantage of DirichletBC(..., method='pointwise') and do something like:

import numpy as np

class Pinpoint(SubDomain):
    TOL = 1e-3
    def __init__(self, coords):
        self.coords = np.array(coords)
        SubDomain.__init__(self)
    def move(self, coords):
        self.coords[:] = np.array(coords)
    def inside(self, x, on_boundary):
        return np.linalg.norm(x-self.coords) < TOL

pinpoint = Pinpoint([0.,0.,0.])

V = VectorFunctionSpace(mesh, "Lagrange", 1)
bc = DirichletBC(V.sub(0), 0.0, pinpoint, 'pointwise')

Note that this has only meaning for V being continuous piecewise linear functions. Constructor Pinpoint.__init__ takes an coords argument (of type tuple, list or np.array) allowing of specification of coordinates where you want to apply BC. Also method Pinpoint.move can be useful to change this point without a need of creating new DirichletBC instance - especially when mesh moves.

Another possibility is to directly fiddle with Tensors' entries.

$\endgroup$
  • $\begingroup$ I will try this function but am little bit clear about your message. thnks $\endgroup$ – nickrocks May 22 '13 at 11:54
  • $\begingroup$ Hi Jan..I tried to work with your idea but little bit confused..I have 88 vertices, 190 cells (3D) in my mesh file and suppose i want to fix x displacement of vertices 84 can you just assist me to do that ?? $\endgroup$ – nickrocks May 23 '13 at 11:12
  • $\begingroup$ @nickrocks: You need to provide more information. I guess you want to apply DirichletBC on some Function - what FunctionSpace is it from? How are these Dirichlet vertices defined? $\endgroup$ – Jan Blechta May 23 '13 at 11:21
  • $\begingroup$ I clarify little bit more..I imported a mesh file which contains vertices (88) and cells in order.. mesh = Mesh("thick_plate_tet_1stOrder.xml") # Define function space and basis functions V = VectorFunctionSpace(mesh, "Lagrange", 1) U = TrialFunction(V) v = TestFunction(V) and I want to fix x direction of displacement only for a particular node (0,0,0) defined in meshfile....how should I do that?? $\endgroup$ – nickrocks May 23 '13 at 11:32
  • $\begingroup$ @nickrocks: Suggestion added to my answer. $\endgroup$ – Jan Blechta May 23 '13 at 12:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.