dirichletBC function giving inconsistent results (in my hands)
Dear all,
I'm running dolfin 1.0.0+ (from dolfin.__version__) and found that the marking of a dirichlet boundary condition for a geometry I am studying (http://
I'm sure I'm doing something foolhardy here, so if you have any quick advice, it would be very appreciated
Thanks!
pete
-------
from dolfin import *
import numpy as np
class General(SubDomain):
def inside(
return on_boundary
class InsideBoundary(
def inside(self, x, on_boundary):
dist = np.linalg.
state = (dist < (1.25*self.radius) and on_boundary)
if(on_boundary and state==1):
print "IB: %f %f %f %f %f %d %d" % (dist, self.radius, x[0],x[
self.ctr+=state
return state
# simple sphere
#mesh = UnitSphere(4)
#c = mesh.coordinates() * 2.25 - [9.5,0,0]
#mesh.coordinat
# my mesh
mesh = Mesh("temp_
# subdom:
subdomains = MeshFunction(
# assign all boundaries with one marker (will be overwritten later)
g = General()
g.mark(
# assign 'left' sphere with a different marker
insideBoundary1 = InsideBoundary()
insideBoundary1
insideBoundary1
insideBoundary1
insideBoundary1
print "Number of times vertices marked %d" % insideBoundary1.ctr
# apply bc cond
V = FunctionSpace(
bc1 = DirichletBC(V, Constant(1), subdomains,1)
bc2 = DirichletBC(V, Constant(0), subdomains,2)
# also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
# also same behaviorbc2 = DirichletBC(V, Constant(0), General())
# view bc2
marked2 = Function(V)
viewbc= DirichletBC(V, Constant(1), subdomains,2)
viewbc.
whereIdx = np.where(
n = np.size(whereIdx)
print "Number of vertices satisfying bc2 %d " % n
#print mesh.coordinates()[ (whereIdx[0])[0],:]
File("test.pvd") << marked2
# solve PDE to verify F'd up
u = TrialFunction(V)
v = TestFunction(V)
# for simple sphere
# form = inner(grad(
# for my geom
form = inner(grad(
A = lhs(form)
L = rhs(form)
x = Function(V)
solve(A=
#print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
File("out.pvd") << x
Question information
- Language:
- English Edit question
- Status:
- Answered
- For:
- DOLFIN Edit question
- Assignee:
- No assignee Edit question
- Last query:
- Last reply:
|
Revision history for this message
|
#1 |
To start, I'd suggest that you write the MeshFunction 'subdomains' to a .pvd file to visualise the marked facets.
|
Revision history for this message
|
#2 |
Thanks Garth,
Sorry for the simple question, but what's the best way to print the
subdomains object?
I tried
File("x.pvd") << subdomains
but the resulting pvd file does not appear correctly in paraview (may have
lost information pertaining to where the surface is located; also the
values range from 0-2e9, not 1-2 as written in the code).
Thanks so far!
pete
On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
<email address hidden>> wrote:
> Your question #214679 on DOLFIN changed:
> https:/
>
> Status: Open => Answered
>
> Garth Wells proposed the following answer:
> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
> .pvd file to visualise the marked facets.
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https:/
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https:/
>
> You received this question notification because you asked the question.
>
--
=======
Peter Kekenes-Huskey, Ph.D.
Postdoctoral Researcher
UCSD Dept of Pharmacology
4205 Urey Hall
9500 Gilman Avenue
La Jolla, CA 92093
http://
=======
|
Revision history for this message
|
#3 |
On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
<email address hidden> wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Status: Answered => Open
>
> Peter Kekenes-Huskey is still having a problem:
> Thanks Garth,
> Sorry for the simple question, but what's the best way to print the
> subdomains object?
> I tried
> File("x.pvd") << subdomains
> but the resulting pvd file does not appear correctly in paraview (may have
> lost information pertaining to where the surface is located; also the
> values range from 0-2e9, not 1-2 as written in the code).
>
The above should work. Which version of DOLFIN are you using?
Garth
> Thanks so far!
> pete
>
>
> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> <email address hidden>> wrote:
>
>> Your question #214679 on DOLFIN changed:
>> https:/
>>
>> Status: Open => Answered
>>
>> Garth Wells proposed the following answer:
>> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
>> .pvd file to visualise the marked facets.
>>
>> --
>> If this answers your question, please go to the following page to let us
>> know that it is solved:
>> https:/
>>
>> If you still need help, you can reply to this email or go to the
>> following page to enter your feedback:
>> https:/
>>
>> You received this question notification because you asked the question.
>>
>
>
> --
> =======
> Peter Kekenes-Huskey, Ph.D.
> Postdoctoral Researcher
> UCSD Dept of Pharmacology
> 4205 Urey Hall
> 9500 Gilman Avenue
> La Jolla, CA 92093
> http://
> =======
>
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
|
Revision history for this message
|
#4 |
The large value is the "no value marker", use a threshold filter in
paraview to hide it.
Martin
Den 20. nov. 2012 18:21 skrev "Garth Wells" <
<email address hidden>> følgende:
> Question #214679 on DOLFIN changed:
> https:/
>
> Status: Open => Answered
>
> Garth Wells proposed the following answer:
> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
> <email address hidden> wrote:
> > Question #214679 on DOLFIN changed:
> > https:/
> >
> > Status: Answered => Open
> >
> > Peter Kekenes-Huskey is still having a problem:
> > Thanks Garth,
> > Sorry for the simple question, but what's the best way to print the
> > subdomains object?
> > I tried
> > File("x.pvd") << subdomains
> > but the resulting pvd file does not appear correctly in paraview (may
> have
> > lost information pertaining to where the surface is located; also the
> > values range from 0-2e9, not 1-2 as written in the code).
> >
>
> The above should work. Which version of DOLFIN are you using?
>
> Garth
>
> > Thanks so far!
> > pete
> >
> >
> > On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> > <email address hidden>> wrote:
> >
> >> Your question #214679 on DOLFIN changed:
> >> https:/
> >>
> >> Status: Open => Answered
> >>
> >> Garth Wells proposed the following answer:
> >> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
> >> .pvd file to visualise the marked facets.
> >>
> >> --
> >> If this answers your question, please go to the following page to let us
> >> know that it is solved:
> >>
> https:/
> >>
> >> If you still need help, you can reply to this email or go to the
> >> following page to enter your feedback:
> >> https:/
> >>
> >> You received this question notification because you asked the question.
> >>
> >
> >
> > --
> > =======
> > Peter Kekenes-Huskey, Ph.D.
> > Postdoctoral Researcher
> > UCSD Dept of Pharmacology
> > 4205 Urey Hall
> > 9500 Gilman Avenue
> > La Jolla, CA 92093
> > http://
> > =======
> >
> > You received this question notification because you are a member of
> > DOLFIN Team, which is an answer contact for DOLFIN.
>
> --
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
>
|
Revision history for this message
|
#5 |
Just a quick note, this can also be visualised with the recently
merged plotting (in trunk). A decent result can be had with f.x.
plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
Wireframe is necessary to see beyond the surface (and can also be
toggled with the 'w' key). Inactive (no value) will be
semi-transparent (light grey). Note that the "hide_above" parameter is
undocumented, just because I didn't know there was a legitimate use
for it before now ;)
On 21 November 2012 08:26, Martin Sandve Alnæs
<email address hidden> wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Martin Sandve Alnæs proposed the following answer:
> The large value is the "no value marker", use a threshold filter in
> paraview to hide it.
>
> Martin
> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
> <email address hidden>> følgende:
>
>> Question #214679 on DOLFIN changed:
>> https:/
>>
>> Status: Open => Answered
>>
>> Garth Wells proposed the following answer:
>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>> <email address hidden> wrote:
>> > Question #214679 on DOLFIN changed:
>> > https:/
>> >
>> > Status: Answered => Open
>> >
>> > Peter Kekenes-Huskey is still having a problem:
>> > Thanks Garth,
>> > Sorry for the simple question, but what's the best way to print the
>> > subdomains object?
>> > I tried
>> > File("x.pvd") << subdomains
>> > but the resulting pvd file does not appear correctly in paraview (may
>> have
>> > lost information pertaining to where the surface is located; also the
>> > values range from 0-2e9, not 1-2 as written in the code).
>> >
>>
>> The above should work. Which version of DOLFIN are you using?
>>
>> Garth
>>
>> > Thanks so far!
>> > pete
>> >
>> >
>> > On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>> > <email address hidden>> wrote:
>> >
>> >> Your question #214679 on DOLFIN changed:
>> >> https:/
>> >>
>> >> Status: Open => Answered
>> >>
>> >> Garth Wells proposed the following answer:
>> >> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
>> >> .pvd file to visualise the marked facets.
>> >>
>> >> --
>> >> If this answers your question, please go to the following page to let us
>> >> know that it is solved:
>> >>
>> https:/
>> >>
>> >> If you still need help, you can reply to this email or go to the
>> >> following page to enter your feedback:
>> >> https:/
>> >>
>> >> You received this question notification because you asked the question.
>> >>
>> >
>> >
>> > --
>> > =======
>> > Peter Kekenes-Huskey, Ph.D.
>> > Postdoctoral Researcher
>> > UCSD Dept of Pharmacology
>> > 4205 Urey Hall
>> > 9500 Gilman Avenue
>> > La Jolla, CA 92093
>> > http://
>> > =======
>> >
>> > You received this question notification because you are a member of
>> > DOLFIN Team, which is an answer contact for DOLFIN.
>>
>> --
>> You received this question notification because you are a member of
>> DOLFIN Team, which is an answer contact for DOLFIN.
>>
>
> --
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
|
Revision history for this message
|
#6 |
The hide_above parameter seems very useful. Could you add it to the
docs?
--
Anders
On Wed, Nov 21, 2012 at 08:41:26AM -0000, Joachim Haga wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Joachim Haga proposed the following answer:
> Just a quick note, this can also be visualised with the recently
> merged plotting (in trunk). A decent result can be had with f.x.
>
> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>
> Wireframe is necessary to see beyond the surface (and can also be
> toggled with the 'w' key). Inactive (no value) will be
> semi-transparent (light grey). Note that the "hide_above" parameter is
> undocumented, just because I didn't know there was a legitimate use
> for it before now ;)
>
>
>
> On 21 November 2012 08:26, Martin Sandve Alnæs
> <email address hidden> wrote:
> > Question #214679 on DOLFIN changed:
> > https:/
> >
> > Martin Sandve Alnæs proposed the following answer:
> > The large value is the "no value marker", use a threshold filter in
> > paraview to hide it.
> >
> > Martin
> > Den 20. nov. 2012 18:21 skrev "Garth Wells" <
> > <email address hidden>> følgende:
> >
> >> Question #214679 on DOLFIN changed:
> >> https:/
> >>
> >> Status: Open => Answered
> >>
> >> Garth Wells proposed the following answer:
> >> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
> >> <email address hidden> wrote:
> >> > Question #214679 on DOLFIN changed:
> >> > https:/
> >> >
> >> > Status: Answered => Open
> >> >
> >> > Peter Kekenes-Huskey is still having a problem:
> >> > Thanks Garth,
> >> > Sorry for the simple question, but what's the best way to print the
> >> > subdomains object?
> >> > I tried
> >> > File("x.pvd") << subdomains
> >> > but the resulting pvd file does not appear correctly in paraview (may
> >> have
> >> > lost information pertaining to where the surface is located; also the
> >> > values range from 0-2e9, not 1-2 as written in the code).
> >> >
> >>
> >> The above should work. Which version of DOLFIN are you using?
> >>
> >> Garth
> >>
> >> > Thanks so far!
> >> > pete
> >> >
> >> >
> >> > On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> >> > <email address hidden>> wrote:
> >> >
> >> >> Your question #214679 on DOLFIN changed:
> >> >> https:/
> >> >>
> >> >> Status: Open => Answered
> >> >>
> >> >> Garth Wells proposed the following answer:
> >> >> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
> >> >> .pvd file to visualise the marked facets.
> >> >>
> >> >>
> >> https:/
> >> >>
> >> >> If you still need help, you can reply to this email or go to the
> >> >> following page to enter your feedback:
> >> >> https:/
> >> >>
> >> >> You received this question notification because you asked the question.
> >> >>
> >> >
> >> >
> >> >
> >> > You received this question notification because you are a member of
> >> > DOLFIN Team, which is an answer contact for DOLFIN.
> >>
> >>
> >
>
|
Revision history for this message
|
#7 |
Very cool joachim. Not sure Pete has access to the latest DOLFIN though :)
Also I cannot reproduce the error you get. With the slightly modified
script below, the correct sphere get marked, and the solution of the
Laplacian look correct.
Johan
from dolfin import *
import numpy as np
class General(SubDomain):
def inside(
return on_boundary
class InsideBoundary(
def inside(self, x, on_boundary):
dist = np.linalg.
state = (dist < (1.25*self.radius) and on_boundary)
if(on_boundary and state==1):
print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
x[0],x[
self.ctr+=state
return state
# my mesh
mesh = Mesh("temp_
# subdom:
# Latest dolfin use sizet instead of uint
#subdomains = FacetFunction(
subdomains = FacetFunction(
# assign all boundaries with one marker (will be overwritten later)
g = General()
g.mark(subdomains, 1)
print "Num faces on boundary:", np.sum(
# assign 'left' sphere with a different marker
insideBoundary1 = InsideBoundary()
insideBoundary1
insideBoundary1
insideBoundary1
insideBoundary1
print "Number of times vertices marked %d" % insideBoundary1.ctr
# apply bc cond
V = FunctionSpace(mesh, "CG", 1)
bc1 = DirichletBC(V, Constant(0), subdomains, 1)
bc2 = DirichletBC(V, Constant(1), subdomains, 2)
# also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
# also same behaviorbc2 = DirichletBC(V, Constant(0), General())
# view bc1
marked1 = Function(V)
marked1.vector()[:] = 0.0
viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
viewbc1.
# view bc2
marked2 = Function(V)
marked2.vector()[:] = 0.0
viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
viewbc2.
whereIdx = np.where(
n = np.size(whereIdx)
print "Number of vertices satisfying bc1 %d " % n
whereIdx = np.where(
n = np.size(whereIdx)
print "Number of vertices satisfying bc2 %d " % n
File("test.pvd") << marked2
# solve PDE to verify F'd up
u = TrialFunction(V)
v = TestFunction(V)
# for my geom
form = inner(grad(u), grad(v))*dx(1)
A = lhs(form)
L = rhs(form)
x = Function(V)
solve(A==L, x, [bc1,bc2])
#print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
#print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
File("out.pvd") << x
#plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
On 11/21/2012 09:41 AM, Joachim Haga wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Joachim Haga proposed the following answer:
> Just a quick note, this can also be visualised with the recently
> merged plotting (in trunk). A decent result can be had with f.x.
>
> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>
> Wireframe is necessary to see beyond the surface (and can also be
> toggled with the 'w' key). Inactive (no value) will be
> semi-transparent (light grey). Note that the "hide_above" parameter is
> undocumented, just because I didn't know there was a legitimate use
> for it before now ;)
>
>
>
> On 21 November 2012 08:26, Martin Sandve Alnæs
> <email address hidden> wrote:
>> Question #214679 on DOLFIN changed:
>> https:/
>>
>> Martin Sandve Alnæs proposed the following answer:
>> The large value is the "no value marker", use a threshold filter in
>> paraview to hide it.
>>
>> Martin
>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
>> <email address hidden>> følgende:
>>
>>> Question #214679 on DOLFIN changed:
>>> https:/
>>>
>>> Status: Open => Answered
>>>
>>> Garth Wells proposed the following answer:
>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>>> <email address hidden> wrote:
>>>> Question #214679 on DOLFIN changed:
>>>> https:/
>>>>
>>>> Status: Answered => Open
>>>>
>>>> Peter Kekenes-Huskey is still having a problem:
>>>> Thanks Garth,
>>>> Sorry for the simple question, but what's the best way to print the
>>>> subdomains object?
>>>> I tried
>>>> File("x.pvd") << subdomains
>>>> but the resulting pvd file does not appear correctly in paraview (may
>>> have
>>>> lost information pertaining to where the surface is located; also the
>>>> values range from 0-2e9, not 1-2 as written in the code).
>>>>
>>>
>>> The above should work. Which version of DOLFIN are you using?
>>>
>>> Garth
>>>
>>>> Thanks so far!
>>>> pete
>>>>
>>>>
>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>>>> <email address hidden>> wrote:
>>>>
>>>>> Your question #214679 on DOLFIN changed:
>>>>> https:/
>>>>>
>>>>> Status: Open => Answered
>>>>>
>>>>> Garth Wells proposed the following answer:
>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
>>>>> .pvd file to visualise the marked facets.
>>>>>
>>>>> --
>>>>> If this answers your question, please go to the following page to let us
>>>>> know that it is solved:
>>>>>
>>> https:/
>>>>>
>>>>> If you still need help, you can reply to this email or go to the
>>>>> following page to enter your feedback:
>>>>> https:/
>>>>>
>>>>> You received this question notification because you asked the question.
>>>>>
>>>>
>>>>
>>>> --
>>>> =======
>>>> Peter Kekenes-Huskey, Ph.D.
>>>> Postdoctoral Researcher
>>>> UCSD Dept of Pharmacology
>>>> 4205 Urey Hall
>>>> 9500 Gilman Avenue
>>>> La Jolla, CA 92093
>>>> http://
>>>> =======
>>>>
>>>> You received this question notification because you are a member of
>>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>
>>> --
>>> You received this question notification because you are a member of
>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>
>>
>> --
>> You received this question notification because you are a member of
>> DOLFIN Team, which is an answer contact for DOLFIN.
>
|
Revision history for this message
|
#8 |
Will do!
j.
On 21 November 2012 11:16, Anders Logg
<email address hidden> wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Anders Logg proposed the following answer:
> The hide_above parameter seems very useful. Could you add it to the
> docs?
>
> --
> Anders
>
>
> On Wed, Nov 21, 2012 at 08:41:26AM -0000, Joachim Haga wrote:
>> Question #214679 on DOLFIN changed:
>> https:/
>>
>> Joachim Haga proposed the following answer:
>> Just a quick note, this can also be visualised with the recently
>> merged plotting (in trunk). A decent result can be had with f.x.
>>
>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>>
>> Wireframe is necessary to see beyond the surface (and can also be
>> toggled with the 'w' key). Inactive (no value) will be
>> semi-transparent (light grey). Note that the "hide_above" parameter is
>> undocumented, just because I didn't know there was a legitimate use
>> for it before now ;)
>>
>>
>>
>> On 21 November 2012 08:26, Martin Sandve Alnæs
>> <email address hidden> wrote:
>> > Question #214679 on DOLFIN changed:
>> > https:/
>> >
>> > Martin Sandve Alnæs proposed the following answer:
>> > The large value is the "no value marker", use a threshold filter in
>> > paraview to hide it.
>> >
>> > Martin
>> > Den 20. nov. 2012 18:21 skrev "Garth Wells" <
>> > <email address hidden>> følgende:
>> >
>> >> Question #214679 on DOLFIN changed:
>> >> https:/
>> >>
>> >> Status: Open => Answered
>> >>
>> >> Garth Wells proposed the following answer:
>> >> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>> >> <email address hidden> wrote:
>> >> > Question #214679 on DOLFIN changed:
>> >> > https:/
>> >> >
>> >> > Status: Answered => Open
>> >> >
>> >> > Peter Kekenes-Huskey is still having a problem:
>> >> > Thanks Garth,
>> >> > Sorry for the simple question, but what's the best way to print the
>> >> > subdomains object?
>> >> > I tried
>> >> > File("x.pvd") << subdomains
>> >> > but the resulting pvd file does not appear correctly in paraview (may
>> >> have
>> >> > lost information pertaining to where the surface is located; also the
>> >> > values range from 0-2e9, not 1-2 as written in the code).
>> >> >
>> >>
>> >> The above should work. Which version of DOLFIN are you using?
>> >>
>> >> Garth
>> >>
>> >> > Thanks so far!
>> >> > pete
>> >> >
>> >> >
>> >> > On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>> >> > <email address hidden>> wrote:
>> >> >
>> >> >> Your question #214679 on DOLFIN changed:
>> >> >> https:/
>> >> >>
>> >> >> Status: Open => Answered
>> >> >>
>> >> >> Garth Wells proposed the following answer:
>> >> >> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
>> >> >> .pvd file to visualise the marked facets.
>> >> >>
>> >> >>
>> >> https:/
>> >> >>
>> >> >> If you still need help, you can reply to this email or go to the
>> >> >> following page to enter your feedback:
>> >> >> https:/
>> >> >>
>> >> >> You received this question notification because you asked the question.
>> >> >>
>> >> >
>> >> >
>> >> >
>> >> > You received this question notification because you are a member of
>> >> > DOLFIN Team, which is an answer contact for DOLFIN.
>> >>
>> >>
>> >
>>
>
> --
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
|
Revision history for this message
|
#9 |
Hi Garth,
I'm using version 1.0 it seems:
>python -c "import dolfin; print dolfin.__version__"
1.0.0
On Tue, Nov 20, 2012 at 9:21 AM, Garth Wells <
<email address hidden>> wrote:
> Your question #214679 on DOLFIN changed:
> https:/
>
> Status: Open => Answered
>
> Garth Wells proposed the following answer:
> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
> <email address hidden> wrote:
> > Question #214679 on DOLFIN changed:
> > https:/
> >
> > Status: Answered => Open
> >
> > Peter Kekenes-Huskey is still having a problem:
> > Thanks Garth,
> > Sorry for the simple question, but what's the best way to print the
> > subdomains object?
> > I tried
> > File("x.pvd") << subdomains
> > but the resulting pvd file does not appear correctly in paraview (may
> have
> > lost information pertaining to where the surface is located; also the
> > values range from 0-2e9, not 1-2 as written in the code).
> >
>
> The above should work. Which version of DOLFIN are you using?
>
> Garth
>
> > Thanks so far!
> > pete
> >
> >
> > On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> > <email address hidden>> wrote:
> >
> >> Your question #214679 on DOLFIN changed:
> >> https:/
> >>
> >> Status: Open => Answered
> >>
> >> Garth Wells proposed the following answer:
> >> To start, I'd suggest that you write the MeshFunction 'subdomains' to a
> >> .pvd file to visualise the marked facets.
> >>
> >> --
> >> If this answers your question, please go to the following page to let us
> >> know that it is solved:
> >>
> https:/
> >>
> >> If you still need help, you can reply to this email or go to the
> >> following page to enter your feedback:
> >> https:/
> >>
> >> You received this question notification because you asked the question.
> >>
> >
> >
> > --
> > =======
> > Peter Kekenes-Huskey, Ph.D.
> > Postdoctoral Researcher
> > UCSD Dept of Pharmacology
> > 4205 Urey Hall
> > 9500 Gilman Avenue
> > La Jolla, CA 92093
> > http://
> > =======
> >
> > You received this question notification because you are a member of
> > DOLFIN Team, which is an answer contact for DOLFIN.
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https:/
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https:/
>
> You received this question notification because you asked the question.
>
--
=======
Peter Kekenes-Huskey, Ph.D.
Postdoctoral Researcher
UCSD Dept of Pharmacology
4205 Urey Hall
9500 Gilman Avenue
La Jolla, CA 92093
http://
=======
|
Revision history for this message
|
#10 |
Thanks! I'll check this out (literally) over the break sometime. We have
Thanksgiving here now
pete
On Wed, Nov 21, 2012 at 12:41 AM, Joachim Haga <
<email address hidden>> wrote:
> Your question #214679 on DOLFIN changed:
> https:/
>
> Joachim Haga proposed the following answer:
> Just a quick note, this can also be visualised with the recently
> merged plotting (in trunk). A decent result can be had with f.x.
>
> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>
> Wireframe is necessary to see beyond the surface (and can also be
> toggled with the 'w' key). Inactive (no value) will be
> semi-transparent (light grey). Note that the "hide_above" parameter is
> undocumented, just because I didn't know there was a legitimate use
> for it before now ;)
>
>
>
> On 21 November 2012 08:26, Martin Sandve Alnæs
> <email address hidden> wrote:
> > Question #214679 on DOLFIN changed:
> > https:/
> >
> > Martin Sandve Alnæs proposed the following answer:
> > The large value is the "no value marker", use a threshold filter in
> > paraview to hide it.
> >
> > Martin
> > Den 20. nov. 2012 18:21 skrev "Garth Wells" <
> > <email address hidden>> følgende:
> >
> >> Question #214679 on DOLFIN changed:
> >> https:/
> >>
> >> Status: Open => Answered
> >>
> >> Garth Wells proposed the following answer:
> >> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
> >> <email address hidden> wrote:
> >> > Question #214679 on DOLFIN changed:
> >> > https:/
> >> >
> >> > Status: Answered => Open
> >> >
> >> > Peter Kekenes-Huskey is still having a problem:
> >> > Thanks Garth,
> >> > Sorry for the simple question, but what's the best way to print the
> >> > subdomains object?
> >> > I tried
> >> > File("x.pvd") << subdomains
> >> > but the resulting pvd file does not appear correctly in paraview (may
> >> have
> >> > lost information pertaining to where the surface is located; also the
> >> > values range from 0-2e9, not 1-2 as written in the code).
> >> >
> >>
> >> The above should work. Which version of DOLFIN are you using?
> >>
> >> Garth
> >>
> >> > Thanks so far!
> >> > pete
> >> >
> >> >
> >> > On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> >> > <email address hidden>> wrote:
> >> >
> >> >> Your question #214679 on DOLFIN changed:
> >> >> https:/
> >> >>
> >> >> Status: Open => Answered
> >> >>
> >> >> Garth Wells proposed the following answer:
> >> >> To start, I'd suggest that you write the MeshFunction 'subdomains'
> to a
> >> >> .pvd file to visualise the marked facets.
> >> >>
> >> >> --
> >> >> If this answers your question, please go to the following page to
> let us
> >> >> know that it is solved:
> >> >>
> >>
> https:/
> >> >>
> >> >> If you still need help, you can reply to this email or go to the
> >> >> following page to enter your feedback:
> >> >> https:/
> >> >>
> >> >> You received this question notification because you asked the
> question.
> >> >>
> >> >
> >> >
> >> > --
> >> > =======
> >> > Peter Kekenes-Huskey, Ph.D.
> >> > Postdoctoral Researcher
> >> > UCSD Dept of Pharmacology
> >> > 4205 Urey Hall
> >> > 9500 Gilman Avenue
> >> > La Jolla, CA 92093
> >> > http://
> >> > =======
> >> >
> >> > You received this question notification because you are a member of
> >> > DOLFIN Team, which is an answer contact for DOLFIN.
> >>
> >> --
> >> You received this question notification because you are a member of
> >> DOLFIN Team, which is an answer contact for DOLFIN.
> >>
> >
> > --
> > You received this question notification because you are a member of
> > DOLFIN Team, which is an answer contact for DOLFIN.
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https:/
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https:/
>
> You received this question notification because you asked the question.
>
--
=======
Peter Kekenes-Huskey, Ph.D.
Postdoctoral Researcher
UCSD Dept of Pharmacology
4205 Urey Hall
9500 Gilman Avenue
La Jolla, CA 92093
http://
=======
|
Revision history for this message
|
#11 |
Johan,
It looks like the 'bug' was that I was using MeshFunction, and not
FacetFunction as you have here. If MeshFunction is used, the number of
marked vertices is like Las Vegas. FacetFunction returns a unique (and
correct) number.
I was under the impression that
MeshFunction(
equivalent - what is the immediate difference (besides a wrong answer in my
case!)
thanks!!
p
On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
<email address hidden>> wrote:
> Your question #214679 on DOLFIN changed:
> https:/
>
> Johan Hake proposed the following answer:
> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
> :)
>
> Also I cannot reproduce the error you get. With the slightly modified
> script below, the correct sphere get marked, and the solution of the
> Laplacian look correct.
>
> Johan
>
> from dolfin import *
> import numpy as np
>
> class General(SubDomain):
> def inside(
> return on_boundary
>
> class InsideBoundary(
> def inside(self, x, on_boundary):
> dist = np.linalg.
> state = (dist < (1.25*self.radius) and on_boundary)
> if(on_boundary and state==1):
> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
> x[0],x[
> self.ctr+=state
> return state
>
> # my mesh
> mesh = Mesh("temp_
>
> # subdom:
> # Latest dolfin use sizet instead of uint
> #subdomains = FacetFunction(
> subdomains = FacetFunction(
>
> # assign all boundaries with one marker (will be overwritten later)
> g = General()
> g.mark(subdomains, 1)
>
> print "Num faces on boundary:", np.sum(
>
> # assign 'left' sphere with a different marker
> insideBoundary1 = InsideBoundary()
> insideBoundary1
> insideBoundary1
> insideBoundary1
> insideBoundary1
>
> print "Number of times vertices marked %d" % insideBoundary1.ctr
>
> # apply bc cond
> V = FunctionSpace(mesh, "CG", 1)
> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
>
> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
>
> # view bc1
> marked1 = Function(V)
> marked1.vector()[:] = 0.0
> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
> viewbc1.
>
> # view bc2
> marked2 = Function(V)
> marked2.vector()[:] = 0.0
> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
> viewbc2.
>
> whereIdx = np.where(
> n = np.size(whereIdx)
> print "Number of vertices satisfying bc1 %d " % n
>
> whereIdx = np.where(
> n = np.size(whereIdx)
> print "Number of vertices satisfying bc2 %d " % n
>
> File("test.pvd") << marked2
>
> # solve PDE to verify F'd up
> u = TrialFunction(V)
> v = TestFunction(V)
>
> # for my geom
> form = inner(grad(u), grad(v))*dx(1)
> A = lhs(form)
> L = rhs(form)
> x = Function(V)
> solve(A==L, x, [bc1,bc2])
>
> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
> File("out.pvd") << x
> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>
>
> On 11/21/2012 09:41 AM, Joachim Haga wrote:
> > Question #214679 on DOLFIN changed:
> > https:/
> >
> > Joachim Haga proposed the following answer:
> > Just a quick note, this can also be visualised with the recently
> > merged plotting (in trunk). A decent result can be had with f.x.
> >
> > plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
> >
> > Wireframe is necessary to see beyond the surface (and can also be
> > toggled with the 'w' key). Inactive (no value) will be
> > semi-transparent (light grey). Note that the "hide_above" parameter is
> > undocumented, just because I didn't know there was a legitimate use
> > for it before now ;)
> >
> >
> >
> > On 21 November 2012 08:26, Martin Sandve Alnæs
> > <email address hidden> wrote:
> >> Question #214679 on DOLFIN changed:
> >> https:/
> >>
> >> Martin Sandve Alnæs proposed the following answer:
> >> The large value is the "no value marker", use a threshold filter in
> >> paraview to hide it.
> >>
> >> Martin
> >> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
> >> <email address hidden>> følgende:
> >>
> >>> Question #214679 on DOLFIN changed:
> >>> https:/
> >>>
> >>> Status: Open => Answered
> >>>
> >>> Garth Wells proposed the following answer:
> >>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
> >>> <email address hidden> wrote:
> >>>> Question #214679 on DOLFIN changed:
> >>>> https:/
> >>>>
> >>>> Status: Answered => Open
> >>>>
> >>>> Peter Kekenes-Huskey is still having a problem:
> >>>> Thanks Garth,
> >>>> Sorry for the simple question, but what's the best way to print the
> >>>> subdomains object?
> >>>> I tried
> >>>> File("x.pvd") << subdomains
> >>>> but the resulting pvd file does not appear correctly in paraview (may
> >>> have
> >>>> lost information pertaining to where the surface is located; also the
> >>>> values range from 0-2e9, not 1-2 as written in the code).
> >>>>
> >>>
> >>> The above should work. Which version of DOLFIN are you using?
> >>>
> >>> Garth
> >>>
> >>>> Thanks so far!
> >>>> pete
> >>>>
> >>>>
> >>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> >>>> <email address hidden>> wrote:
> >>>>
> >>>>> Your question #214679 on DOLFIN changed:
> >>>>> https:/
> >>>>>
> >>>>> Status: Open => Answered
> >>>>>
> >>>>> Garth Wells proposed the following answer:
> >>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
> to a
> >>>>> .pvd file to visualise the marked facets.
> >>>>>
> >>>>> --
> >>>>> If this answers your question, please go to the following page to
> let us
> >>>>> know that it is solved:
> >>>>>
> >>>
> https:/
> >>>>>
> >>>>> If you still need help, you can reply to this email or go to the
> >>>>> following page to enter your feedback:
> >>>>> https:/
> >>>>>
> >>>>> You received this question notification because you asked the
> question.
> >>>>>
> >>>>
> >>>>
> >>>> --
> >>>> =======
> >>>> Peter Kekenes-Huskey, Ph.D.
> >>>> Postdoctoral Researcher
> >>>> UCSD Dept of Pharmacology
> >>>> 4205 Urey Hall
> >>>> 9500 Gilman Avenue
> >>>> La Jolla, CA 92093
> >>>> http://
> >>>> =======
> >>>>
> >>>> You received this question notification because you are a member of
> >>>> DOLFIN Team, which is an answer contact for DOLFIN.
> >>>
> >>> --
> >>> You received this question notification because you are a member of
> >>> DOLFIN Team, which is an answer contact for DOLFIN.
> >>>
> >>
> >> --
> >> You received this question notification because you are a member of
> >> DOLFIN Team, which is an answer contact for DOLFIN.
> >
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https:/
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https:/
>
> You received this question notification because you asked the question.
>
--
=======
Peter Kekenes-Huskey, Ph.D.
Postdoctoral Researcher
UCSD Dept of Pharmacology
4205 Urey Hall
9500 Gilman Avenue
La Jolla, CA 92093
http://
=======
|
Revision history for this message
|
#12 |
On 11/21/2012 04:06 PM, Peter Kekenes-Huskey wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Peter Kekenes-Huskey gave more information on the question:
> Johan,
> It looks like the 'bug' was that I was using MeshFunction, and not
> FacetFunction as you have here. If MeshFunction is used, the number of
> marked vertices is like Las Vegas. FacetFunction returns a unique (and
> correct) number.
Then the problem is uninitialized data...
> I was under the impression that
>
> MeshFunction(
> equivalent - what is the immediate difference (besides a wrong answer in my
> case!)
FacetFunction(
For a 3D mesh.
Johan
> thanks!!
> p
>
>
> On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
> <email address hidden>> wrote:
>
>> Your question #214679 on DOLFIN changed:
>> https:/
>>
>> Johan Hake proposed the following answer:
>> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
>> :)
>>
>> Also I cannot reproduce the error you get. With the slightly modified
>> script below, the correct sphere get marked, and the solution of the
>> Laplacian look correct.
>>
>> Johan
>>
>> from dolfin import *
>> import numpy as np
>>
>> class General(SubDomain):
>> def inside(
>> return on_boundary
>>
>> class InsideBoundary(
>> def inside(self, x, on_boundary):
>> dist = np.linalg.
>> state = (dist < (1.25*self.radius) and on_boundary)
>> if(on_boundary and state==1):
>> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
>> x[0],x[
>> self.ctr+=state
>> return state
>>
>> # my mesh
>> mesh = Mesh("temp_
>>
>> # subdom:
>> # Latest dolfin use sizet instead of uint
>> #subdomains = FacetFunction(
>> subdomains = FacetFunction(
>>
>> # assign all boundaries with one marker (will be overwritten later)
>> g = General()
>> g.mark(subdomains, 1)
>>
>> print "Num faces on boundary:", np.sum(
>>
>> # assign 'left' sphere with a different marker
>> insideBoundary1 = InsideBoundary()
>> insideBoundary1
>> insideBoundary1
>> insideBoundary1
>> insideBoundary1
>>
>> print "Number of times vertices marked %d" % insideBoundary1.ctr
>>
>> # apply bc cond
>> V = FunctionSpace(mesh, "CG", 1)
>> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
>> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
>>
>> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
>> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
>>
>> # view bc1
>> marked1 = Function(V)
>> marked1.vector()[:] = 0.0
>> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
>> viewbc1.
>>
>> # view bc2
>> marked2 = Function(V)
>> marked2.vector()[:] = 0.0
>> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
>> viewbc2.
>>
>> whereIdx = np.where(
>> n = np.size(whereIdx)
>> print "Number of vertices satisfying bc1 %d " % n
>>
>> whereIdx = np.where(
>> n = np.size(whereIdx)
>> print "Number of vertices satisfying bc2 %d " % n
>>
>> File("test.pvd") << marked2
>>
>> # solve PDE to verify F'd up
>> u = TrialFunction(V)
>> v = TestFunction(V)
>>
>> # for my geom
>> form = inner(grad(u), grad(v))*dx(1)
>> A = lhs(form)
>> L = rhs(form)
>> x = Function(V)
>> solve(A==L, x, [bc1,bc2])
>>
>> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
>> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
>> File("out.pvd") << x
>> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>>
>>
>> On 11/21/2012 09:41 AM, Joachim Haga wrote:
>>> Question #214679 on DOLFIN changed:
>>> https:/
>>>
>>> Joachim Haga proposed the following answer:
>>> Just a quick note, this can also be visualised with the recently
>>> merged plotting (in trunk). A decent result can be had with f.x.
>>>
>>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>>>
>>> Wireframe is necessary to see beyond the surface (and can also be
>>> toggled with the 'w' key). Inactive (no value) will be
>>> semi-transparent (light grey). Note that the "hide_above" parameter is
>>> undocumented, just because I didn't know there was a legitimate use
>>> for it before now ;)
>>>
>>>
>>>
>>> On 21 November 2012 08:26, Martin Sandve Alnæs
>>> <email address hidden> wrote:
>>>> Question #214679 on DOLFIN changed:
>>>> https:/
>>>>
>>>> Martin Sandve Alnæs proposed the following answer:
>>>> The large value is the "no value marker", use a threshold filter in
>>>> paraview to hide it.
>>>>
>>>> Martin
>>>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
>>>> <email address hidden>> følgende:
>>>>
>>>>> Question #214679 on DOLFIN changed:
>>>>> https:/
>>>>>
>>>>> Status: Open => Answered
>>>>>
>>>>> Garth Wells proposed the following answer:
>>>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>>>>> <email address hidden> wrote:
>>>>>> Question #214679 on DOLFIN changed:
>>>>>> https:/
>>>>>>
>>>>>> Status: Answered => Open
>>>>>>
>>>>>> Peter Kekenes-Huskey is still having a problem:
>>>>>> Thanks Garth,
>>>>>> Sorry for the simple question, but what's the best way to print the
>>>>>> subdomains object?
>>>>>> I tried
>>>>>> File("x.pvd") << subdomains
>>>>>> but the resulting pvd file does not appear correctly in paraview (may
>>>>> have
>>>>>> lost information pertaining to where the surface is located; also the
>>>>>> values range from 0-2e9, not 1-2 as written in the code).
>>>>>>
>>>>>
>>>>> The above should work. Which version of DOLFIN are you using?
>>>>>
>>>>> Garth
>>>>>
>>>>>> Thanks so far!
>>>>>> pete
>>>>>>
>>>>>>
>>>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>>>>>> <email address hidden>> wrote:
>>>>>>
>>>>>>> Your question #214679 on DOLFIN changed:
>>>>>>> https:/
>>>>>>>
>>>>>>> Status: Open => Answered
>>>>>>>
>>>>>>> Garth Wells proposed the following answer:
>>>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
>> to a
>>>>>>> .pvd file to visualise the marked facets.
>>>>>>>
>>>>>>> --
>>>>>>> If this answers your question, please go to the following page to
>> let us
>>>>>>> know that it is solved:
>>>>>>>
>>>>>
>> https:/
>>>>>>>
>>>>>>> If you still need help, you can reply to this email or go to the
>>>>>>> following page to enter your feedback:
>>>>>>> https:/
>>>>>>>
>>>>>>> You received this question notification because you asked the
>> question.
>>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> =======
>>>>>> Peter Kekenes-Huskey, Ph.D.
>>>>>> Postdoctoral Researcher
>>>>>> UCSD Dept of Pharmacology
>>>>>> 4205 Urey Hall
>>>>>> 9500 Gilman Avenue
>>>>>> La Jolla, CA 92093
>>>>>> http://
>>>>>> =======
>>>>>>
>>>>>> You received this question notification because you are a member of
>>>>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>>>
>>>>> --
>>>>> You received this question notification because you are a member of
>>>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>>>
>>>>
>>>> --
>>>> You received this question notification because you are a member of
>>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>
>>
>> --
>> If this answers your question, please go to the following page to let us
>> know that it is solved:
>> https:/
>>
>> If you still need help, you can reply to this email or go to the
>> following page to enter your feedback:
>> https:/
>>
>> You received this question notification because you asked the question.
>>
>
>
|
Revision history for this message
|
#13 |
OK - thanks for the help everyone and for clarifying my questions
best
pete
On Nov 21, 2012, at 7:36 AM, Johan Hake <email address hidden> wrote:
> Your question #214679 on DOLFIN changed:
> https:/
>
> Status: Open => Answered
>
> Johan Hake proposed the following answer:
> On 11/21/2012 04:06 PM, Peter Kekenes-Huskey wrote:
>> Question #214679 on DOLFIN changed:
>> https:/
>>
>> Peter Kekenes-Huskey gave more information on the question:
>> Johan,
>> It looks like the 'bug' was that I was using MeshFunction, and not
>> FacetFunction as you have here. If MeshFunction is used, the number of
>> marked vertices is like Las Vegas. FacetFunction returns a unique (and
>> correct) number.
>
> Then the problem is uninitialized data...
>
>> I was under the impression that
>>
>> MeshFunction(
>> equivalent - what is the immediate difference (besides a wrong answer in my
>> case!)
>
> FacetFunction(
>
> For a 3D mesh.
>
> Johan
>
>> thanks!!
>> p
>>
>>
>> On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
>> <email address hidden>> wrote:
>>
>>> Your question #214679 on DOLFIN changed:
>>> https:/
>>>
>>> Johan Hake proposed the following answer:
>>> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
>>> :)
>>>
>>> Also I cannot reproduce the error you get. With the slightly modified
>>> script below, the correct sphere get marked, and the solution of the
>>> Laplacian look correct.
>>>
>>> Johan
>>>
>>> from dolfin import *
>>> import numpy as np
>>>
>>> class General(SubDomain):
>>> def inside(
>>> return on_boundary
>>>
>>> class InsideBoundary(
>>> def inside(self, x, on_boundary):
>>> dist = np.linalg.
>>> state = (dist < (1.25*self.radius) and on_boundary)
>>> if(on_boundary and state==1):
>>> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
>>> x[0],x[
>>> self.ctr+=state
>>> return state
>>>
>>> # my mesh
>>> mesh = Mesh("temp_
>>>
>>> # subdom:
>>> # Latest dolfin use sizet instead of uint
>>> #subdomains = FacetFunction(
>>> subdomains = FacetFunction(
>>>
>>> # assign all boundaries with one marker (will be overwritten later)
>>> g = General()
>>> g.mark(subdomains, 1)
>>>
>>> print "Num faces on boundary:", np.sum(
>>>
>>> # assign 'left' sphere with a different marker
>>> insideBoundary1 = InsideBoundary()
>>> insideBoundary1
>>> insideBoundary1
>>> insideBoundary1
>>> insideBoundary1
>>>
>>> print "Number of times vertices marked %d" % insideBoundary1.ctr
>>>
>>> # apply bc cond
>>> V = FunctionSpace(mesh, "CG", 1)
>>> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
>>> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
>>>
>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
>>>
>>> # view bc1
>>> marked1 = Function(V)
>>> marked1.vector()[:] = 0.0
>>> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
>>> viewbc1.
>>>
>>> # view bc2
>>> marked2 = Function(V)
>>> marked2.vector()[:] = 0.0
>>> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
>>> viewbc2.
>>>
>>> whereIdx = np.where(
>>> n = np.size(whereIdx)
>>> print "Number of vertices satisfying bc1 %d " % n
>>>
>>> whereIdx = np.where(
>>> n = np.size(whereIdx)
>>> print "Number of vertices satisfying bc2 %d " % n
>>>
>>> File("test.pvd") << marked2
>>>
>>> # solve PDE to verify F'd up
>>> u = TrialFunction(V)
>>> v = TestFunction(V)
>>>
>>> # for my geom
>>> form = inner(grad(u), grad(v))*dx(1)
>>> A = lhs(form)
>>> L = rhs(form)
>>> x = Function(V)
>>> solve(A==L, x, [bc1,bc2])
>>>
>>> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
>>> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
>>> File("out.pvd") << x
>>> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>>>
>>>
>>> On 11/21/2012 09:41 AM, Joachim Haga wrote:
>>>> Question #214679 on DOLFIN changed:
>>>> https:/
>>>>
>>>> Joachim Haga proposed the following answer:
>>>> Just a quick note, this can also be visualised with the recently
>>>> merged plotting (in trunk). A decent result can be had with f.x.
>>>>
>>>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>>>>
>>>> Wireframe is necessary to see beyond the surface (and can also be
>>>> toggled with the 'w' key). Inactive (no value) will be
>>>> semi-transparent (light grey). Note that the "hide_above" parameter is
>>>> undocumented, just because I didn't know there was a legitimate use
>>>> for it before now ;)
>>>>
>>>>
>>>>
>>>> On 21 November 2012 08:26, Martin Sandve Alnæs
>>>> <email address hidden> wrote:
>>>>> Question #214679 on DOLFIN changed:
>>>>> https:/
>>>>>
>>>>> Martin Sandve Alnæs proposed the following answer:
>>>>> The large value is the "no value marker", use a threshold filter in
>>>>> paraview to hide it.
>>>>>
>>>>> Martin
>>>>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
>>>>> <email address hidden>> følgende:
>>>>>
>>>>>> Question #214679 on DOLFIN changed:
>>>>>> https:/
>>>>>>
>>>>>> Status: Open => Answered
>>>>>>
>>>>>> Garth Wells proposed the following answer:
>>>>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>>>>>> <email address hidden> wrote:
>>>>>>> Question #214679 on DOLFIN changed:
>>>>>>> https:/
>>>>>>>
>>>>>>> Status: Answered => Open
>>>>>>>
>>>>>>> Peter Kekenes-Huskey is still having a problem:
>>>>>>> Thanks Garth,
>>>>>>> Sorry for the simple question, but what's the best way to print the
>>>>>>> subdomains object?
>>>>>>> I tried
>>>>>>> File("x.pvd") << subdomains
>>>>>>> but the resulting pvd file does not appear correctly in paraview (may
>>>>>> have
>>>>>>> lost information pertaining to where the surface is located; also the
>>>>>>> values range from 0-2e9, not 1-2 as written in the code).
>>>>>>>
>>>>>>
>>>>>> The above should work. Which version of DOLFIN are you using?
>>>>>>
>>>>>> Garth
>>>>>>
>>>>>>> Thanks so far!
>>>>>>> pete
>>>>>>>
>>>>>>>
>>>>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>>>>>>> <email address hidden>> wrote:
>>>>>>>
>>>>>>>> Your question #214679 on DOLFIN changed:
>>>>>>>> https:/
>>>>>>>>
>>>>>>>> Status: Open => Answered
>>>>>>>>
>>>>>>>> Garth Wells proposed the following answer:
>>>>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
>>> to a
>>>>>>>> .pvd file to visualise the marked facets.
>>>>>>>>
>>>>>>>> --
>>>>>>>> If this answers your question, please go to the following page to
>>> let us
>>>>>>>> know that it is solved:
>>>>>>>>
>>>>>>
>>> https:/
>>>>>>>>
>>>>>>>> If you still need help, you can reply to this email or go to the
>>>>>>>> following page to enter your feedback:
>>>>>>>> https:/
>>>>>>>>
>>>>>>>> You received this question notification because you asked the
>>> question.
>>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> --
>>>>>>> =======
>>>>>>> Peter Kekenes-Huskey, Ph.D.
>>>>>>> Postdoctoral Researcher
>>>>>>> UCSD Dept of Pharmacology
>>>>>>> 4205 Urey Hall
>>>>>>> 9500 Gilman Avenue
>>>>>>> La Jolla, CA 92093
>>>>>>> http://
>>>>>>> =======
>>>>>>>
>>>>>>> You received this question notification because you are a member of
>>>>>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>>>>
>>>>>> --
>>>>>> You received this question notification because you are a member of
>>>>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>>>>
>>>>>
>>>>> --
>>>>> You received this question notification because you are a member of
>>>>> DOLFIN Team, which is an answer contact for DOLFIN.
>>>>
>>>
>>> --
>>> If this answers your question, please go to the following page to let us
>>> know that it is solved:
>>> https:/
>>>
>>> If you still need help, you can reply to this email or go to the
>>> following page to enter your feedback:
>>> https:/
>>>
>>> You received this question notification because you asked the question.
>>>
>>
>>
>
> --
> If this answers your question, please go to the following page to let us
> know that it is solved:
> https:/
>
> If you still need help, you can reply to this email or go to the
> following page to enter your feedback:
> https:/
>
> You received this question notification because you asked the question.
|
Revision history for this message
|
#14 |
I don't like the MeshFunction constructor that takes an initialization
value. It's confusing as illustrated by this bug report. I would
suggest removing it and instead promoting the use of set_all().
--
Anders
On Wed, Nov 21, 2012 at 04:35:53PM -0000, Peter Kekenes-Huskey wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Status: Answered => Open
>
> Peter Kekenes-Huskey is still having a problem:
> OK - thanks for the help everyone and for clarifying my questions
> best
> pete
>
> On Nov 21, 2012, at 7:36 AM, Johan Hake
> <email address hidden> wrote:
>
> > Your question #214679 on DOLFIN changed:
> > https:/
> >
> > Status: Open => Answered
> >
> > Johan Hake proposed the following answer:
> > On 11/21/2012 04:06 PM, Peter Kekenes-Huskey wrote:
> >> Question #214679 on DOLFIN changed:
> >> https:/
> >>
> >> Peter Kekenes-Huskey gave more information on the question:
> >> Johan,
> >> It looks like the 'bug' was that I was using MeshFunction, and not
> >> FacetFunction as you have here. If MeshFunction is used, the number of
> >> marked vertices is like Las Vegas. FacetFunction returns a unique (and
> >> correct) number.
> >
> > Then the problem is uninitialized data...
> >
> >> I was under the impression that
> >>
> >> MeshFunction(
> >> equivalent - what is the immediate difference (besides a wrong answer in my
> >> case!)
> >
> > FacetFunction(
> >
> > For a 3D mesh.
> >
> > Johan
> >
> >> thanks!!
> >> p
> >>
> >>
> >> On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
> >> <email address hidden>> wrote:
> >>
> >>> Your question #214679 on DOLFIN changed:
> >>> https:/
> >>>
> >>> Johan Hake proposed the following answer:
> >>> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
> >>> :)
> >>>
> >>> Also I cannot reproduce the error you get. With the slightly modified
> >>> script below, the correct sphere get marked, and the solution of the
> >>> Laplacian look correct.
> >>>
> >>> Johan
> >>>
> >>> from dolfin import *
> >>> import numpy as np
> >>>
> >>> class General(SubDomain):
> >>> def inside(
> >>> return on_boundary
> >>>
> >>> class InsideBoundary(
> >>> def inside(self, x, on_boundary):
> >>> dist = np.linalg.
> >>> state = (dist < (1.25*self.radius) and on_boundary)
> >>> if(on_boundary and state==1):
> >>> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
> >>> x[0],x[
> >>> self.ctr+=state
> >>> return state
> >>>
> >>> # my mesh
> >>> mesh = Mesh("temp_
> >>>
> >>> # subdom:
> >>> # Latest dolfin use sizet instead of uint
> >>> #subdomains = FacetFunction(
> >>> subdomains = FacetFunction(
> >>>
> >>> # assign all boundaries with one marker (will be overwritten later)
> >>> g = General()
> >>> g.mark(subdomains, 1)
> >>>
> >>> print "Num faces on boundary:", np.sum(
> >>>
> >>> # assign 'left' sphere with a different marker
> >>> insideBoundary1 = InsideBoundary()
> >>> insideBoundary1
> >>> insideBoundary1
> >>> insideBoundary1
> >>> insideBoundary1
> >>>
> >>> print "Number of times vertices marked %d" % insideBoundary1.ctr
> >>>
> >>> # apply bc cond
> >>> V = FunctionSpace(mesh, "CG", 1)
> >>> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
> >>> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
> >>>
> >>> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
> >>> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
> >>>
> >>> # view bc1
> >>> marked1 = Function(V)
> >>> marked1.vector()[:] = 0.0
> >>> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
> >>> viewbc1.
> >>>
> >>> # view bc2
> >>> marked2 = Function(V)
> >>> marked2.vector()[:] = 0.0
> >>> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
> >>> viewbc2.
> >>>
> >>> whereIdx = np.where(
> >>> n = np.size(whereIdx)
> >>> print "Number of vertices satisfying bc1 %d " % n
> >>>
> >>> whereIdx = np.where(
> >>> n = np.size(whereIdx)
> >>> print "Number of vertices satisfying bc2 %d " % n
> >>>
> >>> File("test.pvd") << marked2
> >>>
> >>> # solve PDE to verify F'd up
> >>> u = TrialFunction(V)
> >>> v = TestFunction(V)
> >>>
> >>> # for my geom
> >>> form = inner(grad(u), grad(v))*dx(1)
> >>> A = lhs(form)
> >>> L = rhs(form)
> >>> x = Function(V)
> >>> solve(A==L, x, [bc1,bc2])
> >>>
> >>> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
> >>> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
> >>> File("out.pvd") << x
> >>> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
> >>>
> >>>
> >>> On 11/21/2012 09:41 AM, Joachim Haga wrote:
> >>>> Question #214679 on DOLFIN changed:
> >>>> https:/
> >>>>
> >>>> Joachim Haga proposed the following answer:
> >>>> Just a quick note, this can also be visualised with the recently
> >>>> merged plotting (in trunk). A decent result can be had with f.x.
> >>>>
> >>>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
> >>>>
> >>>> Wireframe is necessary to see beyond the surface (and can also be
> >>>> toggled with the 'w' key). Inactive (no value) will be
> >>>> semi-transparent (light grey). Note that the "hide_above" parameter is
> >>>> undocumented, just because I didn't know there was a legitimate use
> >>>> for it before now ;)
> >>>>
> >>>>
> >>>>
> >>>> On 21 November 2012 08:26, Martin Sandve Alnæs
> >>>> <email address hidden> wrote:
> >>>>> Question #214679 on DOLFIN changed:
> >>>>> https:/
> >>>>>
> >>>>> Martin Sandve Alnæs proposed the following answer:
> >>>>> The large value is the "no value marker", use a threshold filter in
> >>>>> paraview to hide it.
> >>>>>
> >>>>> Martin
> >>>>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
> >>>>> <email address hidden>> følgende:
> >>>>>
> >>>>>> Question #214679 on DOLFIN changed:
> >>>>>> https:/
> >>>>>>
> >>>>>> Status: Open => Answered
> >>>>>>
> >>>>>> Garth Wells proposed the following answer:
> >>>>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
> >>>>>> <email address hidden> wrote:
> >>>>>>> Question #214679 on DOLFIN changed:
> >>>>>>> https:/
> >>>>>>>
> >>>>>>> Status: Answered => Open
> >>>>>>>
> >>>>>>> Peter Kekenes-Huskey is still having a problem:
> >>>>>>> Thanks Garth,
> >>>>>>> Sorry for the simple question, but what's the best way to print the
> >>>>>>> subdomains object?
> >>>>>>> I tried
> >>>>>>> File("x.pvd") << subdomains
> >>>>>>> but the resulting pvd file does not appear correctly in paraview (may
> >>>>>> have
> >>>>>>> lost information pertaining to where the surface is located; also the
> >>>>>>> values range from 0-2e9, not 1-2 as written in the code).
> >>>>>>>
> >>>>>>
> >>>>>> The above should work. Which version of DOLFIN are you using?
> >>>>>>
> >>>>>> Garth
> >>>>>>
> >>>>>>> Thanks so far!
> >>>>>>> pete
> >>>>>>>
> >>>>>>>
> >>>>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> >>>>>>> <email address hidden>> wrote:
> >>>>>>>
> >>>>>>>> Your question #214679 on DOLFIN changed:
> >>>>>>>> https:/
> >>>>>>>>
> >>>>>>>> Status: Open => Answered
> >>>>>>>>
> >>>>>>>> Garth Wells proposed the following answer:
> >>>>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
> >>> to a
> >>>>>>>> .pvd file to visualise the marked facets.
> >>>>>>>>
> >>> let us
> >>>>>>>> know that it is solved:
> >>>>>>>>
> >>>>>>
> >>> https:/
> >>>>>>>>
> >>>>>>>> If you still need help, you can reply to this email or go to the
> >>>>>>>> following page to enter your feedback:
> >>>>>>>> https:/
> >>>>>>>>
> >>>>>>>> You received this question notification because you asked the
> >>> question.
> >>>>>>>>
> >>>>>>>
> >>>>>>>
> >>>>>>
> >>>>>
> >>>>
> >>>
> >>
> >>
> >
>
|
Revision history for this message
|
#15 |
On Wed, Nov 21, 2012 at 10:15 PM, Anders Logg
<email address hidden> wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Status: Open => Answered
>
> Anders Logg proposed the following answer:
> I don't like the MeshFunction constructor that takes an initialization
> value. It's confusing as illustrated by this bug report. I would
> suggest removing it and instead promoting the use of set_all().
>
I don't think that its confusing and I find it helpful.
Garth
> --
> Anders
>
>
> On Wed, Nov 21, 2012 at 04:35:53PM -0000, Peter Kekenes-Huskey wrote:
>> Question #214679 on DOLFIN changed:
>> https:/
>>
>> Status: Answered => Open
>>
>> Peter Kekenes-Huskey is still having a problem:
>> OK - thanks for the help everyone and for clarifying my questions
>> best
>> pete
>>
>> On Nov 21, 2012, at 7:36 AM, Johan Hake
>> <email address hidden> wrote:
>>
>> > Your question #214679 on DOLFIN changed:
>> > https:/
>> >
>> > Status: Open => Answered
>> >
>> > Johan Hake proposed the following answer:
>> > On 11/21/2012 04:06 PM, Peter Kekenes-Huskey wrote:
>> >> Question #214679 on DOLFIN changed:
>> >> https:/
>> >>
>> >> Peter Kekenes-Huskey gave more information on the question:
>> >> Johan,
>> >> It looks like the 'bug' was that I was using MeshFunction, and not
>> >> FacetFunction as you have here. If MeshFunction is used, the number of
>> >> marked vertices is like Las Vegas. FacetFunction returns a unique (and
>> >> correct) number.
>> >
>> > Then the problem is uninitialized data...
>> >
>> >> I was under the impression that
>> >>
>> >> MeshFunction(
>> >> equivalent - what is the immediate difference (besides a wrong answer in my
>> >> case!)
>> >
>> > FacetFunction(
>> >
>> > For a 3D mesh.
>> >
>> > Johan
>> >
>> >> thanks!!
>> >> p
>> >>
>> >>
>> >> On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
>> >> <email address hidden>> wrote:
>> >>
>> >>> Your question #214679 on DOLFIN changed:
>> >>> https:/
>> >>>
>> >>> Johan Hake proposed the following answer:
>> >>> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
>> >>> :)
>> >>>
>> >>> Also I cannot reproduce the error you get. With the slightly modified
>> >>> script below, the correct sphere get marked, and the solution of the
>> >>> Laplacian look correct.
>> >>>
>> >>> Johan
>> >>>
>> >>> from dolfin import *
>> >>> import numpy as np
>> >>>
>> >>> class General(SubDomain):
>> >>> def inside(
>> >>> return on_boundary
>> >>>
>> >>> class InsideBoundary(
>> >>> def inside(self, x, on_boundary):
>> >>> dist = np.linalg.
>> >>> state = (dist < (1.25*self.radius) and on_boundary)
>> >>> if(on_boundary and state==1):
>> >>> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
>> >>> x[0],x[
>> >>> self.ctr+=state
>> >>> return state
>> >>>
>> >>> # my mesh
>> >>> mesh = Mesh("temp_
>> >>>
>> >>> # subdom:
>> >>> # Latest dolfin use sizet instead of uint
>> >>> #subdomains = FacetFunction(
>> >>> subdomains = FacetFunction(
>> >>>
>> >>> # assign all boundaries with one marker (will be overwritten later)
>> >>> g = General()
>> >>> g.mark(subdomains, 1)
>> >>>
>> >>> print "Num faces on boundary:", np.sum(
>> >>>
>> >>> # assign 'left' sphere with a different marker
>> >>> insideBoundary1 = InsideBoundary()
>> >>> insideBoundary1
>> >>> insideBoundary1
>> >>> insideBoundary1
>> >>> insideBoundary1
>> >>>
>> >>> print "Number of times vertices marked %d" % insideBoundary1.ctr
>> >>>
>> >>> # apply bc cond
>> >>> V = FunctionSpace(mesh, "CG", 1)
>> >>> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
>> >>> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
>> >>>
>> >>> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
>> >>> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
>> >>>
>> >>> # view bc1
>> >>> marked1 = Function(V)
>> >>> marked1.vector()[:] = 0.0
>> >>> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
>> >>> viewbc1.
>> >>>
>> >>> # view bc2
>> >>> marked2 = Function(V)
>> >>> marked2.vector()[:] = 0.0
>> >>> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
>> >>> viewbc2.
>> >>>
>> >>> whereIdx = np.where(
>> >>> n = np.size(whereIdx)
>> >>> print "Number of vertices satisfying bc1 %d " % n
>> >>>
>> >>> whereIdx = np.where(
>> >>> n = np.size(whereIdx)
>> >>> print "Number of vertices satisfying bc2 %d " % n
>> >>>
>> >>> File("test.pvd") << marked2
>> >>>
>> >>> # solve PDE to verify F'd up
>> >>> u = TrialFunction(V)
>> >>> v = TestFunction(V)
>> >>>
>> >>> # for my geom
>> >>> form = inner(grad(u), grad(v))*dx(1)
>> >>> A = lhs(form)
>> >>> L = rhs(form)
>> >>> x = Function(V)
>> >>> solve(A==L, x, [bc1,bc2])
>> >>>
>> >>> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
>> >>> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
>> >>> File("out.pvd") << x
>> >>> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>> >>>
>> >>>
>> >>> On 11/21/2012 09:41 AM, Joachim Haga wrote:
>> >>>> Question #214679 on DOLFIN changed:
>> >>>> https:/
>> >>>>
>> >>>> Joachim Haga proposed the following answer:
>> >>>> Just a quick note, this can also be visualised with the recently
>> >>>> merged plotting (in trunk). A decent result can be had with f.x.
>> >>>>
>> >>>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>> >>>>
>> >>>> Wireframe is necessary to see beyond the surface (and can also be
>> >>>> toggled with the 'w' key). Inactive (no value) will be
>> >>>> semi-transparent (light grey). Note that the "hide_above" parameter is
>> >>>> undocumented, just because I didn't know there was a legitimate use
>> >>>> for it before now ;)
>> >>>>
>> >>>>
>> >>>>
>> >>>> On 21 November 2012 08:26, Martin Sandve Alnæs
>> >>>> <email address hidden> wrote:
>> >>>>> Question #214679 on DOLFIN changed:
>> >>>>> https:/
>> >>>>>
>> >>>>> Martin Sandve Alnæs proposed the following answer:
>> >>>>> The large value is the "no value marker", use a threshold filter in
>> >>>>> paraview to hide it.
>> >>>>>
>> >>>>> Martin
>> >>>>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
>> >>>>> <email address hidden>> følgende:
>> >>>>>
>> >>>>>> Question #214679 on DOLFIN changed:
>> >>>>>> https:/
>> >>>>>>
>> >>>>>> Status: Open => Answered
>> >>>>>>
>> >>>>>> Garth Wells proposed the following answer:
>> >>>>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>> >>>>>> <email address hidden> wrote:
>> >>>>>>> Question #214679 on DOLFIN changed:
>> >>>>>>> https:/
>> >>>>>>>
>> >>>>>>> Status: Answered => Open
>> >>>>>>>
>> >>>>>>> Peter Kekenes-Huskey is still having a problem:
>> >>>>>>> Thanks Garth,
>> >>>>>>> Sorry for the simple question, but what's the best way to print the
>> >>>>>>> subdomains object?
>> >>>>>>> I tried
>> >>>>>>> File("x.pvd") << subdomains
>> >>>>>>> but the resulting pvd file does not appear correctly in paraview (may
>> >>>>>> have
>> >>>>>>> lost information pertaining to where the surface is located; also the
>> >>>>>>> values range from 0-2e9, not 1-2 as written in the code).
>> >>>>>>>
>> >>>>>>
>> >>>>>> The above should work. Which version of DOLFIN are you using?
>> >>>>>>
>> >>>>>> Garth
>> >>>>>>
>> >>>>>>> Thanks so far!
>> >>>>>>> pete
>> >>>>>>>
>> >>>>>>>
>> >>>>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>> >>>>>>> <email address hidden>> wrote:
>> >>>>>>>
>> >>>>>>>> Your question #214679 on DOLFIN changed:
>> >>>>>>>> https:/
>> >>>>>>>>
>> >>>>>>>> Status: Open => Answered
>> >>>>>>>>
>> >>>>>>>> Garth Wells proposed the following answer:
>> >>>>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
>> >>> to a
>> >>>>>>>> .pvd file to visualise the marked facets.
>> >>>>>>>>
>> >>> let us
>> >>>>>>>> know that it is solved:
>> >>>>>>>>
>> >>>>>>
>> >>> https:/
>> >>>>>>>>
>> >>>>>>>> If you still need help, you can reply to this email or go to the
>> >>>>>>>> following page to enter your feedback:
>> >>>>>>>> https:/
>> >>>>>>>>
>> >>>>>>>> You received this question notification because you asked the
>> >>> question.
>> >>>>>>>>
>> >>>>>>>
>> >>>>>>>
>> >>>>>>
>> >>>>>
>> >>>>
>> >>>
>> >>
>> >>
>> >
>>
>
> --
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
|
Revision history for this message
|
#16 |
On 11/21/2012 11:31 PM, Garth Wells wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Garth Wells proposed the following answer:
> On Wed, Nov 21, 2012 at 10:15 PM, Anders Logg
> <email address hidden> wrote:
>> Question #214679 on DOLFIN changed:
>> https:/
>>
>> Status: Open => Answered
>>
>> Anders Logg proposed the following answer:
>> I don't like the MeshFunction constructor that takes an initialization
>> value. It's confusing as illustrated by this bug report. I would
>> suggest removing it and instead promoting the use of set_all().
>>
>
> I don't think that its confusing and I find it helpful.
Me neither, and I also think it is very useful. We could however, always
initialize the MeshFunction to 0. That is what this question is really
about.
Johan
> Garth
>
>> --
>> Anders
>>
>>
>> On Wed, Nov 21, 2012 at 04:35:53PM -0000, Peter Kekenes-Huskey wrote:
>>> Question #214679 on DOLFIN changed:
>>> https:/
>>>
>>> Status: Answered => Open
>>>
>>> Peter Kekenes-Huskey is still having a problem:
>>> OK - thanks for the help everyone and for clarifying my questions
>>> best
>>> pete
>>>
>>> On Nov 21, 2012, at 7:36 AM, Johan Hake
>>> <email address hidden> wrote:
>>>
>>>> Your question #214679 on DOLFIN changed:
>>>> https:/
>>>>
>>>> Status: Open => Answered
>>>>
>>>> Johan Hake proposed the following answer:
>>>> On 11/21/2012 04:06 PM, Peter Kekenes-Huskey wrote:
>>>>> Question #214679 on DOLFIN changed:
>>>>> https:/
>>>>>
>>>>> Peter Kekenes-Huskey gave more information on the question:
>>>>> Johan,
>>>>> It looks like the 'bug' was that I was using MeshFunction, and not
>>>>> FacetFunction as you have here. If MeshFunction is used, the number of
>>>>> marked vertices is like Las Vegas. FacetFunction returns a unique (and
>>>>> correct) number.
>>>>
>>>> Then the problem is uninitialized data...
>>>>
>>>>> I was under the impression that
>>>>>
>>>>> MeshFunction(
>>>>> equivalent - what is the immediate difference (besides a wrong answer in my
>>>>> case!)
>>>>
>>>> FacetFunction(
>>>>
>>>> For a 3D mesh.
>>>>
>>>> Johan
>>>>
>>>>> thanks!!
>>>>> p
>>>>>
>>>>>
>>>>> On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
>>>>> <email address hidden>> wrote:
>>>>>
>>>>>> Your question #214679 on DOLFIN changed:
>>>>>> https:/
>>>>>>
>>>>>> Johan Hake proposed the following answer:
>>>>>> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
>>>>>> :)
>>>>>>
>>>>>> Also I cannot reproduce the error you get. With the slightly modified
>>>>>> script below, the correct sphere get marked, and the solution of the
>>>>>> Laplacian look correct.
>>>>>>
>>>>>> Johan
>>>>>>
>>>>>> from dolfin import *
>>>>>> import numpy as np
>>>>>>
>>>>>> class General(SubDomain):
>>>>>> def inside(
>>>>>> return on_boundary
>>>>>>
>>>>>> class InsideBoundary(
>>>>>> def inside(self, x, on_boundary):
>>>>>> dist = np.linalg.
>>>>>> state = (dist < (1.25*self.radius) and on_boundary)
>>>>>> if(on_boundary and state==1):
>>>>>> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
>>>>>> x[0],x[
>>>>>> self.ctr+=state
>>>>>> return state
>>>>>>
>>>>>> # my mesh
>>>>>> mesh = Mesh("temp_
>>>>>>
>>>>>> # subdom:
>>>>>> # Latest dolfin use sizet instead of uint
>>>>>> #subdomains = FacetFunction(
>>>>>> subdomains = FacetFunction(
>>>>>>
>>>>>> # assign all boundaries with one marker (will be overwritten later)
>>>>>> g = General()
>>>>>> g.mark(subdomains, 1)
>>>>>>
>>>>>> print "Num faces on boundary:", np.sum(
>>>>>>
>>>>>> # assign 'left' sphere with a different marker
>>>>>> insideBoundary1 = InsideBoundary()
>>>>>> insideBoundary1
>>>>>> insideBoundary1
>>>>>> insideBoundary1
>>>>>> insideBoundary1
>>>>>>
>>>>>> print "Number of times vertices marked %d" % insideBoundary1.ctr
>>>>>>
>>>>>> # apply bc cond
>>>>>> V = FunctionSpace(mesh, "CG", 1)
>>>>>> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
>>>>>> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
>>>>>>
>>>>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
>>>>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
>>>>>>
>>>>>> # view bc1
>>>>>> marked1 = Function(V)
>>>>>> marked1.vector()[:] = 0.0
>>>>>> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
>>>>>> viewbc1.
>>>>>>
>>>>>> # view bc2
>>>>>> marked2 = Function(V)
>>>>>> marked2.vector()[:] = 0.0
>>>>>> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
>>>>>> viewbc2.
>>>>>>
>>>>>> whereIdx = np.where(
>>>>>> n = np.size(whereIdx)
>>>>>> print "Number of vertices satisfying bc1 %d " % n
>>>>>>
>>>>>> whereIdx = np.where(
>>>>>> n = np.size(whereIdx)
>>>>>> print "Number of vertices satisfying bc2 %d " % n
>>>>>>
>>>>>> File("test.pvd") << marked2
>>>>>>
>>>>>> # solve PDE to verify F'd up
>>>>>> u = TrialFunction(V)
>>>>>> v = TestFunction(V)
>>>>>>
>>>>>> # for my geom
>>>>>> form = inner(grad(u), grad(v))*dx(1)
>>>>>> A = lhs(form)
>>>>>> L = rhs(form)
>>>>>> x = Function(V)
>>>>>> solve(A==L, x, [bc1,bc2])
>>>>>>
>>>>>> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
>>>>>> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
>>>>>> File("out.pvd") << x
>>>>>> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>>>>>>
>>>>>>
>>>>>> On 11/21/2012 09:41 AM, Joachim Haga wrote:
>>>>>>> Question #214679 on DOLFIN changed:
>>>>>>> https:/
>>>>>>>
>>>>>>> Joachim Haga proposed the following answer:
>>>>>>> Just a quick note, this can also be visualised with the recently
>>>>>>> merged plotting (in trunk). A decent result can be had with f.x.
>>>>>>>
>>>>>>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>>>>>>>
>>>>>>> Wireframe is necessary to see beyond the surface (and can also be
>>>>>>> toggled with the 'w' key). Inactive (no value) will be
>>>>>>> semi-transparent (light grey). Note that the "hide_above" parameter is
>>>>>>> undocumented, just because I didn't know there was a legitimate use
>>>>>>> for it before now ;)
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On 21 November 2012 08:26, Martin Sandve Alnæs
>>>>>>> <email address hidden> wrote:
>>>>>>>> Question #214679 on DOLFIN changed:
>>>>>>>> https:/
>>>>>>>>
>>>>>>>> Martin Sandve Alnæs proposed the following answer:
>>>>>>>> The large value is the "no value marker", use a threshold filter in
>>>>>>>> paraview to hide it.
>>>>>>>>
>>>>>>>> Martin
>>>>>>>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
>>>>>>>> <email address hidden>> følgende:
>>>>>>>>
>>>>>>>>> Question #214679 on DOLFIN changed:
>>>>>>>>> https:/
>>>>>>>>>
>>>>>>>>> Status: Open => Answered
>>>>>>>>>
>>>>>>>>> Garth Wells proposed the following answer:
>>>>>>>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>>>>>>>>> <email address hidden> wrote:
>>>>>>>>>> Question #214679 on DOLFIN changed:
>>>>>>>>>> https:/
>>>>>>>>>>
>>>>>>>>>> Status: Answered => Open
>>>>>>>>>>
>>>>>>>>>> Peter Kekenes-Huskey is still having a problem:
>>>>>>>>>> Thanks Garth,
>>>>>>>>>> Sorry for the simple question, but what's the best way to print the
>>>>>>>>>> subdomains object?
>>>>>>>>>> I tried
>>>>>>>>>> File("x.pvd") << subdomains
>>>>>>>>>> but the resulting pvd file does not appear correctly in paraview (may
>>>>>>>>> have
>>>>>>>>>> lost information pertaining to where the surface is located; also the
>>>>>>>>>> values range from 0-2e9, not 1-2 as written in the code).
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>> The above should work. Which version of DOLFIN are you using?
>>>>>>>>>
>>>>>>>>> Garth
>>>>>>>>>
>>>>>>>>>> Thanks so far!
>>>>>>>>>> pete
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>>>>>>>>>> <email address hidden>> wrote:
>>>>>>>>>>
>>>>>>>>>>> Your question #214679 on DOLFIN changed:
>>>>>>>>>>> https:/
>>>>>>>>>>>
>>>>>>>>>>> Status: Open => Answered
>>>>>>>>>>>
>>>>>>>>>>> Garth Wells proposed the following answer:
>>>>>>>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
>>>>>> to a
>>>>>>>>>>> .pvd file to visualise the marked facets.
>>>>>>>>>>>
>>>>>> let us
>>>>>>>>>>> know that it is solved:
>>>>>>>>>>>
>>>>>>>>>
>>>>>> https:/
>>>>>>>>>>>
>>>>>>>>>>> If you still need help, you can reply to this email or go to the
>>>>>>>>>>> following page to enter your feedback:
>>>>>>>>>>> https:/
>>>>>>>>>>>
>>>>>>>>>>> You received this question notification because you asked the
>>>>>> question.
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>
>>>
>>
>> --
>> You received this question notification because you are a member of
>> DOLFIN Team, which is an answer contact for DOLFIN.
>
|
Revision history for this message
|
#17 |
On Thu, Nov 22, 2012 at 06:45:54AM -0000, Johan Hake wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Johan Hake proposed the following answer:
> On 11/21/2012 11:31 PM, Garth Wells wrote:
> > Question #214679 on DOLFIN changed:
> > https:/
> >
> > Garth Wells proposed the following answer:
> > On Wed, Nov 21, 2012 at 10:15 PM, Anders Logg
> > <email address hidden> wrote:
> >> Question #214679 on DOLFIN changed:
> >> https:/
> >>
> >> Status: Open => Answered
> >>
> >> Anders Logg proposed the following answer:
> >> I don't like the MeshFunction constructor that takes an initialization
> >> value. It's confusing as illustrated by this bug report. I would
> >> suggest removing it and instead promoting the use of set_all().
> >>
> >
> > I don't think that its confusing and I find it helpful.
ok, I won't insist on removing it.
> Me neither, and I also think it is very useful. We could however, always
> initialize the MeshFunction to 0. That is what this question is really
> about.
Are there any objections to or problems with initializing
MeshFunctions with T(0)?
--
Anders
> Johan
>
> > Garth
> >
> >>
> >>
> >> On Wed, Nov 21, 2012 at 04:35:53PM -0000, Peter Kekenes-Huskey wrote:
> >>> Question #214679 on DOLFIN changed:
> >>> https:/
> >>>
> >>> Status: Answered => Open
> >>>
> >>> Peter Kekenes-Huskey is still having a problem:
> >>> OK - thanks for the help everyone and for clarifying my questions
> >>> best
> >>> pete
> >>>
> >>> On Nov 21, 2012, at 7:36 AM, Johan Hake
> >>> <email address hidden> wrote:
> >>>
> >>>> Your question #214679 on DOLFIN changed:
> >>>> https:/
> >>>>
> >>>> Status: Open => Answered
> >>>>
> >>>> Johan Hake proposed the following answer:
> >>>> On 11/21/2012 04:06 PM, Peter Kekenes-Huskey wrote:
> >>>>> Question #214679 on DOLFIN changed:
> >>>>> https:/
> >>>>>
> >>>>> Peter Kekenes-Huskey gave more information on the question:
> >>>>> Johan,
> >>>>> It looks like the 'bug' was that I was using MeshFunction, and not
> >>>>> FacetFunction as you have here. If MeshFunction is used, the number of
> >>>>> marked vertices is like Las Vegas. FacetFunction returns a unique (and
> >>>>> correct) number.
> >>>>
> >>>> Then the problem is uninitialized data...
> >>>>
> >>>>> I was under the impression that
> >>>>>
> >>>>> MeshFunction(
> >>>>> equivalent - what is the immediate difference (besides a wrong answer in my
> >>>>> case!)
> >>>>
> >>>> FacetFunction(
> >>>>
> >>>> For a 3D mesh.
> >>>>
> >>>> Johan
> >>>>
> >>>>> thanks!!
> >>>>> p
> >>>>>
> >>>>>
> >>>>> On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
> >>>>> <email address hidden>> wrote:
> >>>>>
> >>>>>> Your question #214679 on DOLFIN changed:
> >>>>>> https:/
> >>>>>>
> >>>>>> Johan Hake proposed the following answer:
> >>>>>> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
> >>>>>> :)
> >>>>>>
> >>>>>> Also I cannot reproduce the error you get. With the slightly modified
> >>>>>> script below, the correct sphere get marked, and the solution of the
> >>>>>> Laplacian look correct.
> >>>>>>
> >>>>>> Johan
> >>>>>>
> >>>>>> from dolfin import *
> >>>>>> import numpy as np
> >>>>>>
> >>>>>> class General(SubDomain):
> >>>>>> def inside(
> >>>>>> return on_boundary
> >>>>>>
> >>>>>> class InsideBoundary(
> >>>>>> def inside(self, x, on_boundary):
> >>>>>> dist = np.linalg.
> >>>>>> state = (dist < (1.25*self.radius) and on_boundary)
> >>>>>> if(on_boundary and state==1):
> >>>>>> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
> >>>>>> x[0],x[
> >>>>>> self.ctr+=state
> >>>>>> return state
> >>>>>>
> >>>>>> # my mesh
> >>>>>> mesh = Mesh("temp_
> >>>>>>
> >>>>>> # subdom:
> >>>>>> # Latest dolfin use sizet instead of uint
> >>>>>> #subdomains = FacetFunction(
> >>>>>> subdomains = FacetFunction(
> >>>>>>
> >>>>>> # assign all boundaries with one marker (will be overwritten later)
> >>>>>> g = General()
> >>>>>> g.mark(subdomains, 1)
> >>>>>>
> >>>>>> print "Num faces on boundary:", np.sum(
> >>>>>>
> >>>>>> # assign 'left' sphere with a different marker
> >>>>>> insideBoundary1 = InsideBoundary()
> >>>>>> insideBoundary1
> >>>>>> insideBoundary1
> >>>>>> insideBoundary1
> >>>>>> insideBoundary1
> >>>>>>
> >>>>>> print "Number of times vertices marked %d" % insideBoundary1.ctr
> >>>>>>
> >>>>>> # apply bc cond
> >>>>>> V = FunctionSpace(mesh, "CG", 1)
> >>>>>> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
> >>>>>> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
> >>>>>>
> >>>>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
> >>>>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
> >>>>>>
> >>>>>> # view bc1
> >>>>>> marked1 = Function(V)
> >>>>>> marked1.vector()[:] = 0.0
> >>>>>> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
> >>>>>> viewbc1.
> >>>>>>
> >>>>>> # view bc2
> >>>>>> marked2 = Function(V)
> >>>>>> marked2.vector()[:] = 0.0
> >>>>>> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
> >>>>>> viewbc2.
> >>>>>>
> >>>>>> whereIdx = np.where(
> >>>>>> n = np.size(whereIdx)
> >>>>>> print "Number of vertices satisfying bc1 %d " % n
> >>>>>>
> >>>>>> whereIdx = np.where(
> >>>>>> n = np.size(whereIdx)
> >>>>>> print "Number of vertices satisfying bc2 %d " % n
> >>>>>>
> >>>>>> File("test.pvd") << marked2
> >>>>>>
> >>>>>> # solve PDE to verify F'd up
> >>>>>> u = TrialFunction(V)
> >>>>>> v = TestFunction(V)
> >>>>>>
> >>>>>> # for my geom
> >>>>>> form = inner(grad(u), grad(v))*dx(1)
> >>>>>> A = lhs(form)
> >>>>>> L = rhs(form)
> >>>>>> x = Function(V)
> >>>>>> solve(A==L, x, [bc1,bc2])
> >>>>>>
> >>>>>> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
> >>>>>> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
> >>>>>> File("out.pvd") << x
> >>>>>> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
> >>>>>>
> >>>>>>
> >>>>>> On 11/21/2012 09:41 AM, Joachim Haga wrote:
> >>>>>>> Question #214679 on DOLFIN changed:
> >>>>>>> https:/
> >>>>>>>
> >>>>>>> Joachim Haga proposed the following answer:
> >>>>>>> Just a quick note, this can also be visualised with the recently
> >>>>>>> merged plotting (in trunk). A decent result can be had with f.x.
> >>>>>>>
> >>>>>>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
> >>>>>>>
> >>>>>>> Wireframe is necessary to see beyond the surface (and can also be
> >>>>>>> toggled with the 'w' key). Inactive (no value) will be
> >>>>>>> semi-transparent (light grey). Note that the "hide_above" parameter is
> >>>>>>> undocumented, just because I didn't know there was a legitimate use
> >>>>>>> for it before now ;)
> >>>>>>>
> >>>>>>>
> >>>>>>>
> >>>>>>> On 21 November 2012 08:26, Martin Sandve Alnæs
> >>>>>>> <email address hidden> wrote:
> >>>>>>>> Question #214679 on DOLFIN changed:
> >>>>>>>> https:/
> >>>>>>>>
> >>>>>>>> Martin Sandve Alnæs proposed the following answer:
> >>>>>>>> The large value is the "no value marker", use a threshold filter in
> >>>>>>>> paraview to hide it.
> >>>>>>>>
> >>>>>>>> Martin
> >>>>>>>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
> >>>>>>>> <email address hidden>> følgende:
> >>>>>>>>
> >>>>>>>>> Question #214679 on DOLFIN changed:
> >>>>>>>>> https:/
> >>>>>>>>>
> >>>>>>>>> Status: Open => Answered
> >>>>>>>>>
> >>>>>>>>> Garth Wells proposed the following answer:
> >>>>>>>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
> >>>>>>>>> <email address hidden> wrote:
> >>>>>>>>>> Question #214679 on DOLFIN changed:
> >>>>>>>>>> https:/
> >>>>>>>>>>
> >>>>>>>>>> Status: Answered => Open
> >>>>>>>>>>
> >>>>>>>>>> Peter Kekenes-Huskey is still having a problem:
> >>>>>>>>>> Thanks Garth,
> >>>>>>>>>> Sorry for the simple question, but what's the best way to print the
> >>>>>>>>>> subdomains object?
> >>>>>>>>>> I tried
> >>>>>>>>>> File("x.pvd") << subdomains
> >>>>>>>>>> but the resulting pvd file does not appear correctly in paraview (may
> >>>>>>>>> have
> >>>>>>>>>> lost information pertaining to where the surface is located; also the
> >>>>>>>>>> values range from 0-2e9, not 1-2 as written in the code).
> >>>>>>>>>>
> >>>>>>>>>
> >>>>>>>>> The above should work. Which version of DOLFIN are you using?
> >>>>>>>>>
> >>>>>>>>> Garth
> >>>>>>>>>
> >>>>>>>>>> Thanks so far!
> >>>>>>>>>> pete
> >>>>>>>>>>
> >>>>>>>>>>
> >>>>>>>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
> >>>>>>>>>> <email address hidden>> wrote:
> >>>>>>>>>>
> >>>>>>>>>>> Your question #214679 on DOLFIN changed:
> >>>>>>>>>>> https:/
> >>>>>>>>>>>
> >>>>>>>>>>> Status: Open => Answered
> >>>>>>>>>>>
> >>>>>>>>>>> Garth Wells proposed the following answer:
> >>>>>>>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
> >>>>>> to a
> >>>>>>>>>>> .pvd file to visualise the marked facets.
> >>>>>>>>>>>
> >>>>>> let us
> >>>>>>>>>>> know that it is solved:
> >>>>>>>>>>>
> >>>>>>>>>
> >>>>>> https:/
> >>>>>>>>>>>
> >>>>>>>>>>> If you still need help, you can reply to this email or go to the
> >>>>>>>>>>> following page to enter your feedback:
> >>>>>>>>>>> https:/
> >>>>>>>>>>>
> >>>>>>>>>>> You received this question notification because you asked the
> >>>>>> question.
> >>>>>>>>>>>
> >>>>>>>>>>
> >>>>>>>>>>
> >>>>>>>>>
> >>>>>>>>
> >>>>>>>
> >>>>>>
> >>>>>
> >>>>>
> >>>>
> >>>
> >>
> >
>
|
Revision history for this message
|
#18 |
On Thu, Nov 22, 2012 at 6:56 AM, Anders Logg
<email address hidden> wrote:
> Question #214679 on DOLFIN changed:
> https:/
>
> Anders Logg proposed the following answer:
> On Thu, Nov 22, 2012 at 06:45:54AM -0000, Johan Hake wrote:
>> Question #214679 on DOLFIN changed:
>> https:/
>>
>> Johan Hake proposed the following answer:
>> On 11/21/2012 11:31 PM, Garth Wells wrote:
>> > Question #214679 on DOLFIN changed:
>> > https:/
>> >
>> > Garth Wells proposed the following answer:
>> > On Wed, Nov 21, 2012 at 10:15 PM, Anders Logg
>> > <email address hidden> wrote:
>> >> Question #214679 on DOLFIN changed:
>> >> https:/
>> >>
>> >> Status: Open => Answered
>> >>
>> >> Anders Logg proposed the following answer:
>> >> I don't like the MeshFunction constructor that takes an initialization
>> >> value. It's confusing as illustrated by this bug report. I would
>> >> suggest removing it and instead promoting the use of set_all().
>> >>
>> >
>> > I don't think that its confusing and I find it helpful.
>
> ok, I won't insist on removing it.
>
>> Me neither, and I also think it is very useful. We could however, always
>> initialize the MeshFunction to 0. That is what this question is really
>> about.
>
> Are there any objections to or problems with initializing
> MeshFunctions with T(0)?
>
Yes, because MeshFunction can be templated over any object, e.g. a
std::vector of something.
We could use a std::vector for internal storage, which has its own
default initialisation rules. Only problem is using bool because of
some known C++ silliness with std::vectpr<bool>.
Garth
> --
> Anders
>
>
>> Johan
>>
>> > Garth
>> >
>> >>
>> >>
>> >> On Wed, Nov 21, 2012 at 04:35:53PM -0000, Peter Kekenes-Huskey wrote:
>> >>> Question #214679 on DOLFIN changed:
>> >>> https:/
>> >>>
>> >>> Status: Answered => Open
>> >>>
>> >>> Peter Kekenes-Huskey is still having a problem:
>> >>> OK - thanks for the help everyone and for clarifying my questions
>> >>> best
>> >>> pete
>> >>>
>> >>> On Nov 21, 2012, at 7:36 AM, Johan Hake
>> >>> <email address hidden> wrote:
>> >>>
>> >>>> Your question #214679 on DOLFIN changed:
>> >>>> https:/
>> >>>>
>> >>>> Status: Open => Answered
>> >>>>
>> >>>> Johan Hake proposed the following answer:
>> >>>> On 11/21/2012 04:06 PM, Peter Kekenes-Huskey wrote:
>> >>>>> Question #214679 on DOLFIN changed:
>> >>>>> https:/
>> >>>>>
>> >>>>> Peter Kekenes-Huskey gave more information on the question:
>> >>>>> Johan,
>> >>>>> It looks like the 'bug' was that I was using MeshFunction, and not
>> >>>>> FacetFunction as you have here. If MeshFunction is used, the number of
>> >>>>> marked vertices is like Las Vegas. FacetFunction returns a unique (and
>> >>>>> correct) number.
>> >>>>
>> >>>> Then the problem is uninitialized data...
>> >>>>
>> >>>>> I was under the impression that
>> >>>>>
>> >>>>> MeshFunction(
>> >>>>> equivalent - what is the immediate difference (besides a wrong answer in my
>> >>>>> case!)
>> >>>>
>> >>>> FacetFunction(
>> >>>>
>> >>>> For a 3D mesh.
>> >>>>
>> >>>> Johan
>> >>>>
>> >>>>> thanks!!
>> >>>>> p
>> >>>>>
>> >>>>>
>> >>>>> On Wed, Nov 21, 2012 at 3:31 AM, Johan Hake <
>> >>>>> <email address hidden>> wrote:
>> >>>>>
>> >>>>>> Your question #214679 on DOLFIN changed:
>> >>>>>> https:/
>> >>>>>>
>> >>>>>> Johan Hake proposed the following answer:
>> >>>>>> Very cool joachim. Not sure Pete has access to the latest DOLFIN though
>> >>>>>> :)
>> >>>>>>
>> >>>>>> Also I cannot reproduce the error you get. With the slightly modified
>> >>>>>> script below, the correct sphere get marked, and the solution of the
>> >>>>>> Laplacian look correct.
>> >>>>>>
>> >>>>>> Johan
>> >>>>>>
>> >>>>>> from dolfin import *
>> >>>>>> import numpy as np
>> >>>>>>
>> >>>>>> class General(SubDomain):
>> >>>>>> def inside(
>> >>>>>> return on_boundary
>> >>>>>>
>> >>>>>> class InsideBoundary(
>> >>>>>> def inside(self, x, on_boundary):
>> >>>>>> dist = np.linalg.
>> >>>>>> state = (dist < (1.25*self.radius) and on_boundary)
>> >>>>>> if(on_boundary and state==1):
>> >>>>>> print "IB: %f %f %f %f %f %d %d" % (dist, self.radius,
>> >>>>>> x[0],x[
>> >>>>>> self.ctr+=state
>> >>>>>> return state
>> >>>>>>
>> >>>>>> # my mesh
>> >>>>>> mesh = Mesh("temp_
>> >>>>>>
>> >>>>>> # subdom:
>> >>>>>> # Latest dolfin use sizet instead of uint
>> >>>>>> #subdomains = FacetFunction(
>> >>>>>> subdomains = FacetFunction(
>> >>>>>>
>> >>>>>> # assign all boundaries with one marker (will be overwritten later)
>> >>>>>> g = General()
>> >>>>>> g.mark(subdomains, 1)
>> >>>>>>
>> >>>>>> print "Num faces on boundary:", np.sum(
>> >>>>>>
>> >>>>>> # assign 'left' sphere with a different marker
>> >>>>>> insideBoundary1 = InsideBoundary()
>> >>>>>> insideBoundary1
>> >>>>>> insideBoundary1
>> >>>>>> insideBoundary1
>> >>>>>> insideBoundary1
>> >>>>>>
>> >>>>>> print "Number of times vertices marked %d" % insideBoundary1.ctr
>> >>>>>>
>> >>>>>> # apply bc cond
>> >>>>>> V = FunctionSpace(mesh, "CG", 1)
>> >>>>>> bc1 = DirichletBC(V, Constant(0), subdomains, 1)
>> >>>>>> bc2 = DirichletBC(V, Constant(1), subdomains, 2)
>> >>>>>>
>> >>>>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), insideBoundary1)
>> >>>>>> # also same behaviorbc2 = DirichletBC(V, Constant(0), General())
>> >>>>>>
>> >>>>>> # view bc1
>> >>>>>> marked1 = Function(V)
>> >>>>>> marked1.vector()[:] = 0.0
>> >>>>>> viewbc1= DirichletBC(V, Constant(1), subdomains, 1)
>> >>>>>> viewbc1.
>> >>>>>>
>> >>>>>> # view bc2
>> >>>>>> marked2 = Function(V)
>> >>>>>> marked2.vector()[:] = 0.0
>> >>>>>> viewbc2= DirichletBC(V, Constant(1), subdomains, 2)
>> >>>>>> viewbc2.
>> >>>>>>
>> >>>>>> whereIdx = np.where(
>> >>>>>> n = np.size(whereIdx)
>> >>>>>> print "Number of vertices satisfying bc1 %d " % n
>> >>>>>>
>> >>>>>> whereIdx = np.where(
>> >>>>>> n = np.size(whereIdx)
>> >>>>>> print "Number of vertices satisfying bc2 %d " % n
>> >>>>>>
>> >>>>>> File("test.pvd") << marked2
>> >>>>>>
>> >>>>>> # solve PDE to verify F'd up
>> >>>>>> u = TrialFunction(V)
>> >>>>>> v = TestFunction(V)
>> >>>>>>
>> >>>>>> # for my geom
>> >>>>>> form = inner(grad(u), grad(v))*dx(1)
>> >>>>>> A = lhs(form)
>> >>>>>> L = rhs(form)
>> >>>>>> x = Function(V)
>> >>>>>> solve(A==L, x, [bc1,bc2])
>> >>>>>>
>> >>>>>> #print "Value of region near 'absorbing' sphere %f" % x([-9,0,0])
>> >>>>>> #print "Value of region near 'absorbing' sphere %f" % x([-7,0,0])
>> >>>>>> File("out.pvd") << x
>> >>>>>> #plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>> >>>>>>
>> >>>>>>
>> >>>>>> On 11/21/2012 09:41 AM, Joachim Haga wrote:
>> >>>>>>> Question #214679 on DOLFIN changed:
>> >>>>>>> https:/
>> >>>>>>>
>> >>>>>>> Joachim Haga proposed the following answer:
>> >>>>>>> Just a quick note, this can also be visualised with the recently
>> >>>>>>> merged plotting (in trunk). A decent result can be had with f.x.
>> >>>>>>>
>> >>>>>>> plot(subdomains, wireframe=True, hide_above=1e6, interactive=True)
>> >>>>>>>
>> >>>>>>> Wireframe is necessary to see beyond the surface (and can also be
>> >>>>>>> toggled with the 'w' key). Inactive (no value) will be
>> >>>>>>> semi-transparent (light grey). Note that the "hide_above" parameter is
>> >>>>>>> undocumented, just because I didn't know there was a legitimate use
>> >>>>>>> for it before now ;)
>> >>>>>>>
>> >>>>>>>
>> >>>>>>>
>> >>>>>>> On 21 November 2012 08:26, Martin Sandve Alnæs
>> >>>>>>> <email address hidden> wrote:
>> >>>>>>>> Question #214679 on DOLFIN changed:
>> >>>>>>>> https:/
>> >>>>>>>>
>> >>>>>>>> Martin Sandve Alnæs proposed the following answer:
>> >>>>>>>> The large value is the "no value marker", use a threshold filter in
>> >>>>>>>> paraview to hide it.
>> >>>>>>>>
>> >>>>>>>> Martin
>> >>>>>>>> Den 20. nov. 2012 18:21 skrev "Garth Wells" <
>> >>>>>>>> <email address hidden>> følgende:
>> >>>>>>>>
>> >>>>>>>>> Question #214679 on DOLFIN changed:
>> >>>>>>>>> https:/
>> >>>>>>>>>
>> >>>>>>>>> Status: Open => Answered
>> >>>>>>>>>
>> >>>>>>>>> Garth Wells proposed the following answer:
>> >>>>>>>>> On Tue, Nov 20, 2012 at 5:01 PM, Peter Kekenes-Huskey
>> >>>>>>>>> <email address hidden> wrote:
>> >>>>>>>>>> Question #214679 on DOLFIN changed:
>> >>>>>>>>>> https:/
>> >>>>>>>>>>
>> >>>>>>>>>> Status: Answered => Open
>> >>>>>>>>>>
>> >>>>>>>>>> Peter Kekenes-Huskey is still having a problem:
>> >>>>>>>>>> Thanks Garth,
>> >>>>>>>>>> Sorry for the simple question, but what's the best way to print the
>> >>>>>>>>>> subdomains object?
>> >>>>>>>>>> I tried
>> >>>>>>>>>> File("x.pvd") << subdomains
>> >>>>>>>>>> but the resulting pvd file does not appear correctly in paraview (may
>> >>>>>>>>> have
>> >>>>>>>>>> lost information pertaining to where the surface is located; also the
>> >>>>>>>>>> values range from 0-2e9, not 1-2 as written in the code).
>> >>>>>>>>>>
>> >>>>>>>>>
>> >>>>>>>>> The above should work. Which version of DOLFIN are you using?
>> >>>>>>>>>
>> >>>>>>>>> Garth
>> >>>>>>>>>
>> >>>>>>>>>> Thanks so far!
>> >>>>>>>>>> pete
>> >>>>>>>>>>
>> >>>>>>>>>>
>> >>>>>>>>>> On Tue, Nov 20, 2012 at 12:06 AM, Garth Wells <
>> >>>>>>>>>> <email address hidden>> wrote:
>> >>>>>>>>>>
>> >>>>>>>>>>> Your question #214679 on DOLFIN changed:
>> >>>>>>>>>>> https:/
>> >>>>>>>>>>>
>> >>>>>>>>>>> Status: Open => Answered
>> >>>>>>>>>>>
>> >>>>>>>>>>> Garth Wells proposed the following answer:
>> >>>>>>>>>>> To start, I'd suggest that you write the MeshFunction 'subdomains'
>> >>>>>> to a
>> >>>>>>>>>>> .pvd file to visualise the marked facets.
>> >>>>>>>>>>>
>> >>>>>> let us
>> >>>>>>>>>>> know that it is solved:
>> >>>>>>>>>>>
>> >>>>>>>>>
>> >>>>>> https:/
>> >>>>>>>>>>>
>> >>>>>>>>>>> If you still need help, you can reply to this email or go to the
>> >>>>>>>>>>> following page to enter your feedback:
>> >>>>>>>>>>> https:/
>> >>>>>>>>>>>
>> >>>>>>>>>>> You received this question notification because you asked the
>> >>>>>> question.
>> >>>>>>>>>>>
>> >>>>>>>>>>
>> >>>>>>>>>>
>> >>>>>>>>>
>> >>>>>>>>
>> >>>>>>>
>> >>>>>>
>> >>>>>
>> >>>>>
>> >>>>
>> >>>
>> >>
>> >
>>
>
> --
> You received this question notification because you are a member of
> DOLFIN Team, which is an answer contact for DOLFIN.
Can you help with this problem?
Provide an answer of your own, or ask Peter Kekenes-Huskey for more information if necessary.
