[yt-users] Defining a new field, including IF statement

[Britton Smith]

Hi again,

There's a small typo.  This line:
field_data[np.invert(threshold] = 0
is missing the closing bracket:
field_data[np.invert(threshold]] = 0

On Thu, Nov 2, 2017 at 6:00 AM, Britton Smith <brittonsmith at gmail.com>
wrote:

> Hi Yusuke,
>
> This error comes from the fact that your field function is being evaluated
> for an array of elements all at once.  You will likely have portions of the
> array that satisfy both criteria.  You'll need to first define a field
> array and then set the values for the portions above the density threshold
> and then below.  I recommend doing something like the following:
>
> def _SFR_Density(field, data):
>    f_sf = 0.01
>    rho_thresh = 57.5 * Mu * hydrogen_mass
>    # make field array
>    field_data = data["density"].copy()
>    # density threshold
>    threshold = data["density"] > rho_thresh
>    # set values above threshold
>    field_data[threshold] *= f_sf * / (3*pi/32/Grav/data["density"])**(0.5)
>    # set values below
>    field_data[np.invert(threshold] = 0
>    return field_data
>
> The threshold array is a boolean array that can be used to access all the
> values satisfying your criteria all at once.  Then, the np.invert gives you
> all the values not satisfying it.  Hope that helps!
>
> Britton
>
> On Wed, Nov 1, 2017 at 10:46 PM, Yusuke Fujimoto <
> yusuke.fujimoto.jp at gmail.com> wrote:
>
>> Dear yt-users,
>>
>> I am trying to define a new field, SFR_Density, which includes if
>> statements.
>> That is defined by,
>> f_sf * rho / t_ff  if rho > rho_threshold
>> 0                      if rho < rho_threshold
>> (f_sf is star formation efficiency, rho is gas density, and t_ff is free
>> fall time)
>>
>> Here is my definition I tried.
>>
>> def _SFR_Density(field, data):
>>
>>    f_sf = 0.01
>>
>>    rho_thresh = 57.5 * Mu * hydrogen_mass
>>
>>    if data["density"] > rho_thresh:
>>
>>       return f_sf * data["density"] / (3*pi/32/Grav/data["density"])**(
>> 0.5)
>>
>>    else:
>>
>>       return 0.0 * f_sf * data["density"] / (3*pi/32/Grav/data["density"
>> ])**(0.5)
>>
>> add_field("SFR_Density", function=_SFR_Density)
>>
>> But it fails with an error message of "ValueError: The truth value of an
>> array with more than one element is ambiguous. Use a.any() or a.all()"
>> I am sure this is because data["density"] is a list, so I need to specify
>> the exact index number of elements at the if statement. I have tried many
>> tests, but now I am in stuck.
>>
>> Any help would be very appreciated.
>>
>> Thanks,
>> Yusuke
>>
>> _______________________________________________
>> yt-users mailing list
>> yt-users at lists.spacepope.org
>> http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.spacepope.org/pipermail/yt-users-spacepope.org/attachments/20171102/e2474779/attachment-0001.html>