Hi guys, I've been exploring the idea of changing radius to be correct for periodic boxes. Right now it is incorrect; each component (x,y,z) should not have a distance greater than 0.5 * domain_size. The easiest way to do this would be: rx = min( abs( x  center_x ), abs( x  center_x  domain_x) ) I think the best way to do this is to set up a NumPy ufunc http://docs.scipy.org/doc/numpy/user/cinfo.beyondbasics.html#creatingane... that accepts three arrays, along with the center and the domain width, and then returns the radius. What do you all think? Alternatively, I'm thinking maybe just a regular function that gets the arrays and returns one would be better; the ufunc machinery is a bit complicated and I might get confused. Once I come up with it, will somebody be able to look over my work? Matt
Hi guys, Can somebody review this? I chatted with Jeff, and I believe it is correct for periodic boxes, but I need somebody to look at it before I commit it, as it is a rather crucial aspect of yt. This changes the way radius is calculated and should now work with periodic boundary conditions. Here's the new field definition: def _NewRadius(field, data): center = data.get_field_parameter("center") DW = data.pf["DomainRightEdge"]  data.pf["DomainLeftEdge"] radius = na.zeros(data["x"].shape, dtype='float64') for i, ax in enumerate('xyz'): r = na.abs(data[ax]  center[i]) radius += na.minimum(r, na.abs(DW[i]r))**2.0 na.sqrt(radius, radius) return radius the final argument in na.sqrt is the 'out' argument, to avoid temporary arrays. I've committed a sample script here: http://paste.enzotools.org/show/101/ which loads up a data dump, calculates the radius from [0.1, 0.1, 0.1] at all points in the root grid, then plots the output. It gives correct min & max values for the radius. If I can get somebody else to sign off on this, I will commit the necessary changes to _Radius and _ParticleRadius, as well as grid & point selection for spheres. We'll then have a fullyperiodic profiler, which will close #165. Matt On Thu, Apr 16, 2009 at 11:39 AM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
I've been exploring the idea of changing radius to be correct for periodic boxes. Right now it is incorrect; each component (x,y,z) should not have a distance greater than 0.5 * domain_size. The easiest way to do this would be:
rx = min( abs( x  center_x ), abs( x  center_x  domain_x) )
I think the best way to do this is to set up a NumPy ufunc
http://docs.scipy.org/doc/numpy/user/cinfo.beyondbasics.html#creatingane...
that accepts three arrays, along with the center and the domain width, and then returns the radius. What do you all think? Alternatively, I'm thinking maybe just a regular function that gets the arrays and returns one would be better; the ufunc machinery is a bit complicated and I might get confused. Once I come up with it, will somebody be able to look over my work?
Matt
Oops, just as a note, there's an error in the paste  the second plot of "NewRadius" should be of "Radius". Thanks guys... Matt On Fri, Apr 17, 2009 at 7:16 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
Can somebody review this? I chatted with Jeff, and I believe it is correct for periodic boxes, but I need somebody to look at it before I commit it, as it is a rather crucial aspect of yt. This changes the way radius is calculated and should now work with periodic boundary conditions.
Here's the new field definition:
def _NewRadius(field, data): center = data.get_field_parameter("center") DW = data.pf["DomainRightEdge"]  data.pf["DomainLeftEdge"] radius = na.zeros(data["x"].shape, dtype='float64') for i, ax in enumerate('xyz'): r = na.abs(data[ax]  center[i]) radius += na.minimum(r, na.abs(DW[i]r))**2.0 na.sqrt(radius, radius) return radius
the final argument in na.sqrt is the 'out' argument, to avoid temporary arrays.
I've committed a sample script here:
http://paste.enzotools.org/show/101/
which loads up a data dump, calculates the radius from [0.1, 0.1, 0.1] at all points in the root grid, then plots the output. It gives correct min & max values for the radius.
If I can get somebody else to sign off on this, I will commit the necessary changes to _Radius and _ParticleRadius, as well as grid & point selection for spheres. We'll then have a fullyperiodic profiler, which will close #165.
Matt
On Thu, Apr 16, 2009 at 11:39 AM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
I've been exploring the idea of changing radius to be correct for periodic boxes. Right now it is incorrect; each component (x,y,z) should not have a distance greater than 0.5 * domain_size. The easiest way to do this would be:
rx = min( abs( x  center_x ), abs( x  center_x  domain_x) )
I think the best way to do this is to set up a NumPy ufunc
http://docs.scipy.org/doc/numpy/user/cinfo.beyondbasics.html#creatingane...
that accepts three arrays, along with the center and the domain width, and then returns the radius. What do you all think? Alternatively, I'm thinking maybe just a regular function that gets the arrays and returns one would be better; the ufunc machinery is a bit complicated and I might get confused. Once I come up with it, will somebody be able to look over my work?
Matt
I just tried your script with some data I have here. I've attached the images. I can't really tell if this is right or not. I think another way to test this would be to do a projection with a periodic sphere as the source object. If you throw a sphere that overlaps the domain boundary to a projection, it should be quite clear if it's doing what it's supposed to. Does that sound right? Britton On Fri, Apr 17, 2009 at 8:17 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Oops, just as a note, there's an error in the paste  the second plot of "NewRadius" should be of "Radius".
Thanks guys...
Matt
On Fri, Apr 17, 2009 at 7:16 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
Can somebody review this? I chatted with Jeff, and I believe it is correct for periodic boxes, but I need somebody to look at it before I commit it, as it is a rather crucial aspect of yt. This changes the way radius is calculated and should now work with periodic boundary conditions.
Here's the new field definition:
def _NewRadius(field, data): center = data.get_field_parameter("center") DW = data.pf["DomainRightEdge"]  data.pf["DomainLeftEdge"] radius = na.zeros(data["x"].shape, dtype='float64') for i, ax in enumerate('xyz'): r = na.abs(data[ax]  center[i]) radius += na.minimum(r, na.abs(DW[i]r))**2.0 na.sqrt(radius, radius) return radius
the final argument in na.sqrt is the 'out' argument, to avoid temporary arrays.
I've committed a sample script here:
http://paste.enzotools.org/show/101/
which loads up a data dump, calculates the radius from [0.1, 0.1, 0.1] at all points in the root grid, then plots the output. It gives correct min & max values for the radius.
If I can get somebody else to sign off on this, I will commit the necessary changes to _Radius and _ParticleRadius, as well as grid & point selection for spheres. We'll then have a fullyperiodic profiler, which will close #165.
Matt
On Thu, Apr 16, 2009 at 11:39 AM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
I've been exploring the idea of changing radius to be correct for periodic boxes. Right now it is incorrect; each component (x,y,z) should not have a distance greater than 0.5 * domain_size. The easiest way to do this would be:
rx = min( abs( x  center_x ), abs( x  center_x  domain_x) )
I think the best way to do this is to set up a NumPy ufunc
http://docs.scipy.org/doc/numpy/user/cinfo.beyondbasics.html#creatingane...
that accepts three arrays, along with the center and the domain width, and then returns the radius. What do you all think? Alternatively, I'm thinking maybe just a regular function that gets the arrays and returns one would be better; the ufunc machinery is a bit complicated and I might get confused. Once I come up with it, will somebody be able to look over my work?
Matt
_______________________________________________ Ytdev mailing list Ytdev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/ytdevspacepope.org
Hm, this looks wrong to me. Clearly something I have not considered  over here myi mages look like the attached. What's the value for the domain width? This is confusing to me, as I really thought it would work! Can you try with a projection? It should work, right through the domain, and it should tell us what we're looking for. Or, jsut straightup slices of NewRadius. I'll try some stuff too. On Fri, Apr 17, 2009 at 8:48 PM, Britton Smith <brittonsmith@gmail.com> wrote:
I just tried your script with some data I have here. I've attached the images. I can't really tell if this is right or not. I think another way to test this would be to do a projection with a periodic sphere as the source object. If you throw a sphere that overlaps the domain boundary to a projection, it should be quite clear if it's doing what it's supposed to. Does that sound right?
Britton
On Fri, Apr 17, 2009 at 8:17 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Oops, just as a note, there's an error in the paste  the second plot of "NewRadius" should be of "Radius".
Thanks guys...
Matt
On Fri, Apr 17, 2009 at 7:16 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
Can somebody review this? I chatted with Jeff, and I believe it is correct for periodic boxes, but I need somebody to look at it before I commit it, as it is a rather crucial aspect of yt. This changes the way radius is calculated and should now work with periodic boundary conditions.
Here's the new field definition:
def _NewRadius(field, data): center = data.get_field_parameter("center") DW = data.pf["DomainRightEdge"]  data.pf["DomainLeftEdge"] radius = na.zeros(data["x"].shape, dtype='float64') for i, ax in enumerate('xyz'): r = na.abs(data[ax]  center[i]) radius += na.minimum(r, na.abs(DW[i]r))**2.0 na.sqrt(radius, radius) return radius
the final argument in na.sqrt is the 'out' argument, to avoid temporary arrays.
I've committed a sample script here:
http://paste.enzotools.org/show/101/
which loads up a data dump, calculates the radius from [0.1, 0.1, 0.1] at all points in the root grid, then plots the output. It gives correct min & max values for the radius.
If I can get somebody else to sign off on this, I will commit the necessary changes to _Radius and _ParticleRadius, as well as grid & point selection for spheres. We'll then have a fullyperiodic profiler, which will close #165.
Matt
On Thu, Apr 16, 2009 at 11:39 AM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
I've been exploring the idea of changing radius to be correct for periodic boxes. Right now it is incorrect; each component (x,y,z) should not have a distance greater than 0.5 * domain_size. The easiest way to do this would be:
rx = min( abs( x  center_x ), abs( x  center_x  domain_x) )
I think the best way to do this is to set up a NumPy ufunc
http://docs.scipy.org/doc/numpy/user/cinfo.beyondbasics.html#creatingane...
that accepts three arrays, along with the center and the domain width, and then returns the radius. What do you all think? Alternatively, I'm thinking maybe just a regular function that gets the arrays and returns one would be better; the ufunc machinery is a bit complicated and I might get confused. Once I come up with it, will somebody be able to look over my work?
Matt
_______________________________________________ Ytdev mailing list Ytdev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/ytdevspacepope.org
_______________________________________________ Ytdev mailing list Ytdev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/ytdevspacepope.org
Okay, so suddenly it occurs to me  if your domain width is 1.0, then this is in fact just a subset of your domain, and it could be correct. I forgot about parallel root grid IO. What happens if you plot a slice through the center of the domain of the new field? pc = PlotCollection(pf, center=[0.5,0.5,0.5]) pc.add_slice("NewRadius", 0) pc.add_slice("Radius", 0) pc.save("new") That should show us if it's actually giving the right results or not. Apologies; I completely forgot about decomposed root grids when I wrote that script. On Fri, Apr 17, 2009 at 9:01 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hm, this looks wrong to me. Clearly something I have not considered  over here myi mages look like the attached. What's the value for the domain width? This is confusing to me, as I really thought it would work!
Can you try with a projection? It should work, right through the domain, and it should tell us what we're looking for. Or, jsut straightup slices of NewRadius. I'll try some stuff too.
On Fri, Apr 17, 2009 at 8:48 PM, Britton Smith <brittonsmith@gmail.com> wrote:
I just tried your script with some data I have here. I've attached the images. I can't really tell if this is right or not. I think another way to test this would be to do a projection with a periodic sphere as the source object. If you throw a sphere that overlaps the domain boundary to a projection, it should be quite clear if it's doing what it's supposed to. Does that sound right?
Britton
On Fri, Apr 17, 2009 at 8:17 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Oops, just as a note, there's an error in the paste  the second plot of "NewRadius" should be of "Radius".
Thanks guys...
Matt
On Fri, Apr 17, 2009 at 7:16 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
Can somebody review this? I chatted with Jeff, and I believe it is correct for periodic boxes, but I need somebody to look at it before I commit it, as it is a rather crucial aspect of yt. This changes the way radius is calculated and should now work with periodic boundary conditions.
Here's the new field definition:
def _NewRadius(field, data): center = data.get_field_parameter("center") DW = data.pf["DomainRightEdge"]  data.pf["DomainLeftEdge"] radius = na.zeros(data["x"].shape, dtype='float64') for i, ax in enumerate('xyz'): r = na.abs(data[ax]  center[i]) radius += na.minimum(r, na.abs(DW[i]r))**2.0 na.sqrt(radius, radius) return radius
the final argument in na.sqrt is the 'out' argument, to avoid temporary arrays.
I've committed a sample script here:
http://paste.enzotools.org/show/101/
which loads up a data dump, calculates the radius from [0.1, 0.1, 0.1] at all points in the root grid, then plots the output. It gives correct min & max values for the radius.
If I can get somebody else to sign off on this, I will commit the necessary changes to _Radius and _ParticleRadius, as well as grid & point selection for spheres. We'll then have a fullyperiodic profiler, which will close #165.
Matt
On Thu, Apr 16, 2009 at 11:39 AM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
I've been exploring the idea of changing radius to be correct for periodic boxes. Right now it is incorrect; each component (x,y,z) should not have a distance greater than 0.5 * domain_size. The easiest way to do this would be:
rx = min( abs( x  center_x ), abs( x  center_x  domain_x) )
I think the best way to do this is to set up a NumPy ufunc
http://docs.scipy.org/doc/numpy/user/cinfo.beyondbasics.html#creatingane...
that accepts three arrays, along with the center and the domain width, and then returns the radius. What do you all think? Alternatively, I'm thinking maybe just a regular function that gets the arrays and returns one would be better; the ufunc machinery is a bit complicated and I might get confused. Once I come up with it, will somebody be able to look over my work?
Matt
_______________________________________________ Ytdev mailing list Ytdev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/ytdevspacepope.org
_______________________________________________ Ytdev mailing list Ytdev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/ytdevspacepope.org
Hi guys, I did a bunch of testing on this today and satisfied myself it was correct. I also fixed a regression I introduced, clarified the setting of field parameters, and then made AMRSpheres periodic. This was committed in r1267. Code review is welcome! Matt On Fri, Apr 17, 2009 at 9:10 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Okay, so suddenly it occurs to me  if your domain width is 1.0, then this is in fact just a subset of your domain, and it could be correct. I forgot about parallel root grid IO. What happens if you plot a slice through the center of the domain of the new field?
pc = PlotCollection(pf, center=[0.5,0.5,0.5]) pc.add_slice("NewRadius", 0) pc.add_slice("Radius", 0) pc.save("new")
That should show us if it's actually giving the right results or not. Apologies; I completely forgot about decomposed root grids when I wrote that script.
On Fri, Apr 17, 2009 at 9:01 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hm, this looks wrong to me. Clearly something I have not considered  over here myi mages look like the attached. What's the value for the domain width? This is confusing to me, as I really thought it would work!
Can you try with a projection? It should work, right through the domain, and it should tell us what we're looking for. Or, jsut straightup slices of NewRadius. I'll try some stuff too.
On Fri, Apr 17, 2009 at 8:48 PM, Britton Smith <brittonsmith@gmail.com> wrote:
I just tried your script with some data I have here. I've attached the images. I can't really tell if this is right or not. I think another way to test this would be to do a projection with a periodic sphere as the source object. If you throw a sphere that overlaps the domain boundary to a projection, it should be quite clear if it's doing what it's supposed to. Does that sound right?
Britton
On Fri, Apr 17, 2009 at 8:17 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Oops, just as a note, there's an error in the paste  the second plot of "NewRadius" should be of "Radius".
Thanks guys...
Matt
On Fri, Apr 17, 2009 at 7:16 PM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
Can somebody review this? I chatted with Jeff, and I believe it is correct for periodic boxes, but I need somebody to look at it before I commit it, as it is a rather crucial aspect of yt. This changes the way radius is calculated and should now work with periodic boundary conditions.
Here's the new field definition:
def _NewRadius(field, data): center = data.get_field_parameter("center") DW = data.pf["DomainRightEdge"]  data.pf["DomainLeftEdge"] radius = na.zeros(data["x"].shape, dtype='float64') for i, ax in enumerate('xyz'): r = na.abs(data[ax]  center[i]) radius += na.minimum(r, na.abs(DW[i]r))**2.0 na.sqrt(radius, radius) return radius
the final argument in na.sqrt is the 'out' argument, to avoid temporary arrays.
I've committed a sample script here:
http://paste.enzotools.org/show/101/
which loads up a data dump, calculates the radius from [0.1, 0.1, 0.1] at all points in the root grid, then plots the output. It gives correct min & max values for the radius.
If I can get somebody else to sign off on this, I will commit the necessary changes to _Radius and _ParticleRadius, as well as grid & point selection for spheres. We'll then have a fullyperiodic profiler, which will close #165.
Matt
On Thu, Apr 16, 2009 at 11:39 AM, Matthew Turk <matthewturk@gmail.com> wrote:
Hi guys,
I've been exploring the idea of changing radius to be correct for periodic boxes. Right now it is incorrect; each component (x,y,z) should not have a distance greater than 0.5 * domain_size. The easiest way to do this would be:
rx = min( abs( x  center_x ), abs( x  center_x  domain_x) )
I think the best way to do this is to set up a NumPy ufunc
http://docs.scipy.org/doc/numpy/user/cinfo.beyondbasics.html#creatingane...
that accepts three arrays, along with the center and the domain width, and then returns the radius. What do you all think? Alternatively, I'm thinking maybe just a regular function that gets the arrays and returns one would be better; the ufunc machinery is a bit complicated and I might get confused. Once I come up with it, will somebody be able to look over my work?
Matt
_______________________________________________ Ytdev mailing list Ytdev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/ytdevspacepope.org
_______________________________________________ Ytdev mailing list Ytdev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/ytdevspacepope.org
participants (2)

Britton Smith

Matthew Turk