It does not work because it is not guaranteed that there is any facet at x[0]==42.0. You can try finding approximating facets somehow or fit your mesh to what you need. The latter option should be possible with mshr but continue asking such questions on fenicsproject.org/qa, please.
Jan On Thu, 11 Dec 2014 18:31:51 +0100 simone <[email protected]> wrote: > I've tried running this code following your instructions: > > from dolfin import * > from mshr import * > > geom=Rectangle(dolfin.Point(0.0,0.0),dolfin.Point(50.0,50.0)) > mesh=generate_mesh(geom,100,"cgal") > > class Axis42(SubDomain): > def inside(self, x, on_boundary): > return near(x[0], 42.0) > > facets = FacetFunction('size_t', mesh) > facets.set_all(0) > Axis42().mark(facets, 2) > > for i in range(len(facets.array())): > if facets.array()[i]==2: > print ("marked") > > > the compiler returns no print; if I change line 9 in > > return near(x[0],0.0) > > the mark method works, and facets.array() is not equal to 0 for all i. > Really thank you for your help. > Regards > SC > > > Il 11/12/2014 14:58, Jan Blechta ha scritto: > > On Wed, 10 Dec 2014 20:24:08 +0100 > > simone <[email protected]> wrote: > > > >> Hi, > >> thank's for your support; i've tried using your code, that runs > >> without any error, but it still doesn't do what i need. The problem > >> is that with this definition of facets, the mark method of > >> LeftBoundary and RightBoundary objects (as defined in my first > >> mail) doesn't mark my domain. It is able to mark only parts of > >> domain boundary, not subdomain boundary (which is within the > >> domain, not on the boundary ), so if I try to integrate over > >> subdomain boundary it returns me 0.0 because there's no facet over > >> which to integrate. I suppose this problem is due to my > >> implementation error, could you please explain me where i'm wrong? > >> Thank you > > [please, keep [email protected] in CC] > > > > 1. ds = measure on exterior facets > > dS = measure on interior facets > > > > 2. on_boundary in your SubDomains causes that it is only on > > boundary. Example > > > > class Axis42(SubDomain): > > def inside(self, x, on_boundary): > > return near(x[0], 42.0) > > > > defines axis x=42.0 disregarding boundary. > > > > Jan > > > >> Regards > >> SC > >> > >> Il 10/12/2014 15:11, Jan Blechta ha scritto: > >>> Hi, > >>> > >>> we cannot find a problem why non-positive progress bar is to be > >>> created and raising runtime error without a complete code to try. > >>> > >>> But note that you are probably misusing SubMesh. Having solved > >>> some problem on mesh and trying to integrate its solution around > >>> some cells or facets is normally done without SubMesh. Like > >>> > >>> u = Function(V) # V space on mesh > >>> > >>> solve(..., u, ...) > >>> > >>> facets = FacetFunction('size_t', mesh) > >>> LeftBoundary().mark(facets, 2) > >>> RightBoundary().mark(facets, 1) > >>> > >>> # current syntax > >>> F = u*ds(subdomain_id=1, subdomain_data=facets) > >>> F += u*ds(subdomain_id=2, subdomain_data=facets) > >>> > >>> # deprecated synatx > >>> F = u*ds[facets](1) > >>> F += u*ds[facets](2) > >>> > >>> assemble(F) > >>> > >>> i.e. everything is done on the orginal mesh and integration is > >>> performed over on subset of original mesh using cell/facet > >>> markers. No new mesh is ibtroduced using SubMesh, which just > >>> complicates the things. > >>> > >>> Jan > >>> > >>> > >>> On Wed, 10 Dec 2014 02:19:13 +0100 > >>> simone <[email protected]> wrote: > >>> > >>>> Good evening, > >>>> i'm sorry if I rewrite to you but i cannot find a solution in the > >>>> forum neither in the documentation. I have a rectangular domain > >>>> and a smaller rectangular subdomain. After the resolution of a > >>>> differential problem over the domain, i need to evaluate the > >>>> integral of the obtained solution over the left and right sides > >>>> of the subdomain boundary. Called p the solution, this is part > >>>> of my code: > >>>> > >>>> Gamma=FacetFunction("size_t",SubMesh(my_mesh,my_subdomain)) > >>>> > >>>> class LeftBoundary(SubDomain): > >>>> def inside(self, x, on_boundary): > >>>> return on_boundary and near(x[0],x1) #x1=x-coordinate > >>>> of the left side of the subdomain > >>>> > >>>> class RightBoundary(SubDomain): > >>>> def inside(self, x, on_boundary): > >>>> return on_boundary and near(x[0],x2) #x2=x-coordinate > >>>> of the right side of the subdomain > >>>> > >>>> Gamma.set_all(0) > >>>> LeftBoundary().mark(Gamma,2) > >>>> RightBoundary().mark(Gamma,1) > >>>> > >>>> d_gamma=Measure("ds",domain=SubMesh(my_mesh,my_subdomain))[Gamma] > >>>> > >>>> d_p=p*d_gamma(1)-p*d_gamma(0) > >>>> > >>>> direction=assemble(d_p,exterior_facet_domains=Gamma) > >>>> > >>>> > >>>> > >>>> > >>>> If I try to run this code i get this error: > >>>> File "/usr/lib/python2.7/dist-packages/dolfin/cpp/mesh.py", line > >>>> 4461, in mark > >>>> self._mark(*args) > >>>> 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 *** > >>>> *** [email protected] > >>>> *** > >>>> *** Remember to include the error message listed below and, if > >>>> possible, *** include a *minimal* running example to reproduce > >>>> the error. *** > >>>> *** > >>>> ------------------------------------------------------------------------- > >>>> *** Error: Unable to create progress bar. > >>>> *** Reason: Number of steps for progress bar must be positive. > >>>> *** Where: This error was encountered inside Progress.cpp. > >>>> *** Process: unknown > >>>> *** > >>>> *** DOLFIN version: 1.4.0 > >>>> *** Git changeset: unknown > >>>> > >>>> I cannot understand the meaning of the error message, so i don't > >>>> know how to modify the code. > >>>> I would like also to know what happens if the mesh vertexes don't > >>>> match the subdomain boundary, with the subdomain sides > >>>> intersecting mesh cells. I hope I've explained myself, I'm sorry > >>>> for my english Thank you for your support > >>>> Best regards > >>>> Simone Carriero > >>>> _______________________________________________ > >>>> fenics-support mailing list > >>>> [email protected] > >>>> http://fenicsproject.org/mailman/listinfo/fenics-support > > _______________________________________________ > fenics-support mailing list > [email protected] > http://fenicsproject.org/mailman/listinfo/fenics-support _______________________________________________ fenics-support mailing list [email protected] http://fenicsproject.org/mailman/listinfo/fenics-support
