Dear all,
I am hoping to calculate the flux of mass ("mdot") and angular momentum ("Ldot") through a series of concentric spherical surfaces in a simulation snapshot (centered around an absorbing "sink" within the grid).
So far, I've attempted this along these lines:
pf = load(fn) # load data rad = 0.1 sp = pf.h.sphere([0,0,0],1.1*rad) surf = pf.h.surface(sp,"Radius",rad) flux = surf.calculate_flux("x-velocity","y-velocity","z-velocity","Density")
where flux gets assigned to the mass accretion rate in this case. I believe I could do the angular momentum by calculating the flux of a derived variable that is defined as data["SpecificAngularMomentumZ"] * data["Density"] such that the resulting units are angular momentum per volume.
I was wondering if there is a better way to go about this calculation, it seems that others must have tried to look at similar questions. In particular, defining an "isoradius" surface seems to accomplish what I was hoping to do, but there must be a cleaner way.
Thanks so much for your advice,
Morgan
Hi Morgan,
Sorry for the delay in replying. I find this type of functionality really interesting, so I'm quite happy to see you applying it.
On Thu, Oct 24, 2013 at 10:03 AM, Morgan MacLeod morganmacleod@gmail.com wrote:
Dear all,
I am hoping to calculate the flux of mass ("mdot") and angular momentum ("Ldot") through a series of concentric spherical surfaces in a simulation snapshot (centered around an absorbing "sink" within the grid).
So far, I've attempted this along these lines:
pf = load(fn) # load data rad = 0.1 sp = pf.h.sphere([0,0,0],1.1*rad) surf = pf.h.surface(sp,"Radius",rad) flux = surf.calculate_flux("x-velocity","y-velocity","z-velocity","Density")
where flux gets assigned to the mass accretion rate in this case.
To verify, I have re-read the source code and the comments, and I agree with this assessment. The flux here will be computed by calculating the dot product of the local velocity vector with the normal to the identified surface, compute the area of the triangle in every cell, and then multiply rho * area * (normal dot vel). The units here will be scaled to code area units, so you may need to multiply by the conversion factor for cm**2. But, the end result will be g/s.
I believe I could do the angular momentum by calculating the flux of a derived variable that is defined as data["SpecificAngularMomentumZ"] * data["Density"] such that the resulting units are angular momentum per volume.
I think this is where it gets somewhat tricky. I'm not entirely certain that this gives precisely what you are looking for, but it might. So what this will do is give the Z component of the angular momentum, dotted by the normal vector, times the area. My reluctance here is how the L components add in quadrature and where the summation of that should occur. If you do the three fluxes independently:
flux_x = surf.calculate_flux("x-velocity","y-velocity","z-velocity","AngularMomentumX") flux_y = surf.calculate_flux("x-velocity","y-velocity","z-velocity","AngularMomentumY") flux_z = surf.calculate_flux("x-velocity","y-velocity","z-velocity","AngularMomentumZ")
That's not the same as compute the total L locally (i.e., sqrt(Lx^2 + Ly^2 + Lz^2)) and computing the flux of that over the entire surface. Depending on the question you're asking, one might have meaning and the other might not.
I was wondering if there is a better way to go about this calculation, it seems that others must have tried to look at similar questions. In particular, defining an "isoradius" surface seems to accomplish what I was hoping to do, but there must be a cleaner way.
In principle, it is a clean way of doing it -- the other way would be to do a radial profile and look at the total enclosed L at the varying radii and differencing it. The surface object will interpolate to the surface itself, whereas the radial profile will only do binning.
-Matt
Thanks so much for your advice,
Morgan
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Dear Matt,
Thank you for your reply and advice. These are some important points about the subtleties of the different methods (ie radial profiles and binning vs. interpolation to a surface) that I hadn't fully appreciated.
I really like that this list is a place that one can look for inspiration and new analysis techniques as well as normal support and solutions to technical problems. I know I've benefitted from the discussion.
In any case, I'm attaching a script that measures the flux of mass and the z-component of angular momentum by interpolation to a series of concentric spherical surfaces (in a time series of simulation snapshots). In my simulation, these surfaces surround an absorbing "sink" in the grid. I am using yt-3.0 here, although I'm not sure whether what I've done is version-specific. Perhaps this is a useful starting point to someone else on the list.
Morgan
On Fri, Oct 25, 2013 at 4:06 PM, Matthew Turk matthewturk@gmail.com wrote:
Hi Morgan,
Sorry for the delay in replying. I find this type of functionality really interesting, so I'm quite happy to see you applying it.
On Thu, Oct 24, 2013 at 10:03 AM, Morgan MacLeod morganmacleod@gmail.com wrote:
Dear all,
I am hoping to calculate the flux of mass ("mdot") and angular momentum ("Ldot") through a series of concentric spherical surfaces in a
simulation
snapshot (centered around an absorbing "sink" within the grid).
So far, I've attempted this along these lines:
pf = load(fn) # load data rad = 0.1 sp = pf.h.sphere([0,0,0],1.1*rad) surf = pf.h.surface(sp,"Radius",rad) flux =
surf.calculate_flux("x-velocity","y-velocity","z-velocity","Density")
where flux gets assigned to the mass accretion rate in this case.
To verify, I have re-read the source code and the comments, and I agree with this assessment. The flux here will be computed by calculating the dot product of the local velocity vector with the normal to the identified surface, compute the area of the triangle in every cell, and then multiply rho * area * (normal dot vel). The units here will be scaled to code area units, so you may need to multiply by the conversion factor for cm**2. But, the end result will be g/s.
I believe I could do the angular momentum by calculating the flux of a derived variable that is defined as data["SpecificAngularMomentumZ"] * data["Density"] such that the resulting units are angular momentum per volume.
I think this is where it gets somewhat tricky. I'm not entirely certain that this gives precisely what you are looking for, but it might. So what this will do is give the Z component of the angular momentum, dotted by the normal vector, times the area. My reluctance here is how the L components add in quadrature and where the summation of that should occur. If you do the three fluxes independently:
flux_x = surf.calculate_flux("x-velocity","y-velocity","z-velocity","AngularMomentumX") flux_y = surf.calculate_flux("x-velocity","y-velocity","z-velocity","AngularMomentumY") flux_z = surf.calculate_flux("x-velocity","y-velocity","z-velocity","AngularMomentumZ")
That's not the same as compute the total L locally (i.e., sqrt(Lx^2 + Ly^2 + Lz^2)) and computing the flux of that over the entire surface. Depending on the question you're asking, one might have meaning and the other might not.
I was wondering if there is a better way to go about this calculation, it seems that others must have tried to look at similar questions. In particular, defining an "isoradius" surface seems to accomplish what I
was
hoping to do, but there must be a cleaner way.
In principle, it is a clean way of doing it -- the other way would be to do a radial profile and look at the total enclosed L at the varying radii and differencing it. The surface object will interpolate to the surface itself, whereas the radial profile will only do binning.
-Matt
Thanks so much for your advice,
Morgan
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
participants (2)
-
Matthew Turk -
Morgan MacLeod