In [1]:

    

from dolfin import *
import numpy
nn = 2
mesh = UnitCubeMesh(int(nn),int(nn),int(nn))
order = 1
parameters['reorder_dofs_serial'] = False
Magnetic = FunctionSpace(mesh, "N1curl", order)
Lagrange = FunctionSpace(mesh, "CG", order)
print Magnetic.dim()

def boundary(x, on_boundary):
        return on_boundary
W = Magnetic*Lagrange
# bcW = DirichletBC(W.sub(0),Expression(("1.0","1.0","1.0")), boundary)
# bcuW = DirichletBC(W.sub(1), Expression(("1.0")), boundary)


    

    




DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:FFC:Reusing form from cache.
DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:FFC:Reusing form from cache.
DEBUG:UFL:No integrals left after transformation, returning empty form.
DEBUG:FFC:Reusing form from cache.






    




98

In [3]:

    

a = mesh.bounding_box_tree()
B = BoundaryMesh(Magnetic.mesh(),"exterior",False)
B = numpy.sort(B.entity_map(0).array().astype("float_","C"))


    

In [40]:

    

v = Vertex(mesh, 3)


    

In [61]:

    

data = mesh.data()


    

In [63]:

    

data.array("facet_orientation",0)


    

    




---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
<ipython-input-63-05ed6cbbe53d> in <module>()
----> 1 data.array("facet_orientation",0)

RuntimeError: 

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
***     fenics@fenicsproject.org
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to access mesh data.
*** Reason:  Mesh data array named "facet_orientation" does not exist.
*** Where:   This error was encountered inside MeshData.cpp.
*** Process: unknown
*** 
*** DOLFIN version: 1.4.0
*** Git changeset:  
*** -------------------------------------------------------------------------

In [5]:

    

j = 0
for vert in vertices(mesh):
        if vert.index() == B[j]:
            print "boundary"
            j = j+1
        else:
            print "\nvert", vert.index()
            print vert.entities(1)
            for edge in edges(vert):
#                 j = j+1
#                 print j
                print "   edge", edge.index()
Ve


    

    




---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-5-ccbc6f441f30> in <module>()
     11 #                 print j
     12                 print "   edge", edge.index()
---> 13 Ve

NameError: name 'Ve' is not defined





    




boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary

vert 13
[ 0  1  3  6  9 12 15 32 38 53 57 79 82 93]
   edge 0
   edge 1
   edge 3
   edge 6
   edge 9
   edge 12
   edge 15
   edge 32
   edge 38
   edge 53
   edge 57
   edge 79
   edge 82
   edge 93
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary
boundary

In [76]:

    




    

In [80]:

    

E.array().size


    

    
Out[80]:





98

In [83]:

    

E = EdgeFunctionBool(mesh)
VertexBoundary = numpy.zeros(mesh.num_edges())
VertexBoundary[E.array()] = 1


    

In [84]:

    

VertexBoundary


    

    
Out[84]:





array([ 1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.,  1.,  1.,  1.,  1.,  0.,
        0.,  0.,  0.,  1.,  0.,  1.,  1.,  0.,  1.,  1.,  0.,  1.,  0.,
        0.,  1.,  1.,  0.,  1.,  1.,  0.,  1.,  0.,  0.,  1.,  0.,  0.,
        0.,  1.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  1.,  1.,  1.,  1.,
        1.,  1.,  0.,  0.,  1.,  1.,  1.,  1.,  1.,  1.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,  0.,
        0.,  0.,  0.,  0.,  0.,  0.,  0.])

In [59]:

    

SystemAssembler?


    

In [24]:

    

bcW = DirichletBC(Magnetic,Expression(("1.0","1.0","1.0")), boundary)
bcW.get_boundary_values()
bcW.markers()


    

    
Out[24]:





array([  3,   6,   9,  12,  15,  17,  18,  21,  22,  24,  27,  30,  37,
        42,  44,  47,  48,  49,  50,  53,  54,  58,  60,  63,  69,  70,
        72,  75,  77,  79,  80,  83,  85,  86,  88,  91,  99, 101, 103,
       105, 106, 107, 108, 111, 113, 115, 117, 119], dtype=uint64)

In [16]:

    

B.keys()
intlist = [str(x) for x in B.keys()]


    

In [21]:

    

import numpy
x= numpy.array(intlist)
y = x.astype("float_","C")


    

In [22]:

    

y


    

    
Out[22]:





array([  2.,   4.,   5.,   7.,   8.,  10.,  11.,  13.,  14.,  16.,  17.,
        18.,  19.,  20.,  21.,  23.,  24.,  25.,  26.,  27.,  29.,  31.,
        33.,  35.,  37.,  40.,  41.,  42.,  43.,  44.,  45.,  46.,  47.,
        48.,  49.,  51.,  52.,  54.,  55.,  56.,  60.,  61.,  62.,  63.,
        64.,  65.,  67.,  68.,  69.,  70.,  71.,  72.,  74.,  75.,  76.,
        77.,  78.,  80.,  81.,  84.,  85.,  86.,  87.,  88.,  89.,  90.,
        91.,  92.,  94.,  95.,  96.,  97.])