[SciPy-user] += Operator and Slicing of Arrays

[Josh Lawrence]

Emmanuelle,

Thank you for your speedy reply. I apologize if the list gets this and
another message replying to the same thing. Gmail imap seems to be
down so I am having to use webmail. Anyways, I get an error telling me
that for numpy.histogram, the "weights" variable needs to be the same
shape as "a" which I believe corresponds to tri_idx_plus in the code
example you gave me.

Also, the edge_unkown (and hence quantity) has dtype of complex128
(double valued complex). Does this preclude using numpy.histogram?

Thanks!

Josh Lawrence

On Fri, Jun 12, 2009 at 2:28 PM, Emmanuelle
Gouillart<emmanuelle.gouillart at normalesup.org> wrote:
> Hi Josh,
>
> what kind of problem do you have exactly ? Do you have trouble
> implementing the computation you describe, or do you get unexpected
> results?
>
> If I understood well what you want to do, you cannot use directly use +=
> with fancy indexing (quantity[tri_idx_plus,...]) because the repeated
> elements will be incremented just once (see
> http://www.scipy.org/Tentative_NumPy_Tutorial#head-0dffc419afa7d77d51062d40d2d84143db8216c2
> for more details).
>
> However, I think you can solve your problem by using a weighted
> histogram. Using your notations
>
> weights = edge_unknown[:len(tri_idx_plus),...] *
>        basis_p[:len(tri_idx_plus),...]
> histogram_values = np.histogram(tri_idx_plus,
>        np.arange(tri_idx_plus.max() +2), weights=weights)
> unique_plus = np.unique1d(tri_idx_plus)
> quantity[unique_plus,...] = histogram_values[0][unique_plus]
>
> The weighted histogram allows you to make the sums corresponding to each
> triangle.
>
> Here is an example
>
>>>> quantity = np.zeros(10)
>>>> tri_idx_plus = np.array([0, 3, 6, 3, 2, 4, 1, 4])
>>>> weights = 2*tri_idx_plus + 1
>>>> histogram_values = np.histogram(tri_idx_plus,
>        np.arange(tri_idx_plus.max() +2), weights=weights)
>>>> unique_plus = np.unique1d(tri_idx_plus)
>>>> quantity[unique_plus,...] = histogram_values[0][unique_plus]
>
> Actually, this may only work with positive values of weights (not
> checked)...
>
> Please tell us if this meets your needs or not.
>
> Cheers,
>
> Emmanuelle
>
> On Fri, Jun 12, 2009 at 01:14:16PM -0400, Josh Lawrence wrote:
>> Greetings,
>
>> I am in need of some help. I am trying to use the += operator to sum
>> over edge elements on a triangular mesh. Each edge has an unknown
>> associated with it. After solving for the unknowns, I am trying to
>> compute a quantity at the centroid of all triangles in the mesh. On a
>> flat surface, each edge will be connected to either one or two
>> triangles. The orientation of the edge and the normals of the
>> triangles determines whether each triangle attached to the edge is a
>> "plus" or "minus" triangle for that edge. It is possible for one
>> triangle to be referenced three times as a "plus" triangle, three
>> times as a "minus" triangle or any combination of "plus" and
>> "minus" (1 and 2 or 2 and 1, respectively).
>
>> I have a variable tri_idx which relates the edges to the "plus" and
>> "minus" triangles. I then compute the quantity at the centroid for the
>> "plus" triangle and "minus" triangle attached to each edge. An example
>> follows:
>
>> tri_idx_plus = [0 3 6 3 2 4 1 4]
>> tri_idx_minus = [1 2 5 3 6 0 1 4]
>
>> quantity[tri_idx_plus,...] += edge_unknown[:len(tri_idx_plus),...] *
>> basis_p[:len(tri_idx_plus),...]
>> quantity[tri_idx_minus,...] += edge_unknown[:len(tri_idx_minus)...] *
>> basis_m[:len(tri_idx_minus),...]
>
>> where basis_p and basis_m are basis functions that expand the unknown
>> of each edge into a surface function over the "plus" or "minus"
>> triangle.
>
>> I am pretty sure the problem I am encountering is that tri_idx_plus
>> mentions indices 3 and 4 twice and tri_idx_minus contains index 1
>> twice. Is there a way of doing this operation without reverting to
>> looping over each edge (read: not doing this the slow way).
>
>> Thanks in advance!
>
>> Josh Lawrence
>> Ph.D. Student
>> Clemson University
>
>> _______________________________________________
>> SciPy-user mailing list
>> SciPy-user at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-user
> _______________________________________________
> SciPy-user mailing list
> SciPy-user at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>



-- 
Josh Lawrence
Ph.D. Student
Clemson University