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

[Emmanuelle Gouillart]

> It tells me that the weights argument needs to be the same shape as  
> tri_idx_plus. Also, edge_unknown is complex valued. So I have +- real  
> and imaginary components. I do not know much about numpy.histogram...  
> Does the complex dtype preclude it from use?

	Yes, the weight array must have the same shape as the values
array ("a" in the help of np.histogram). This is the case in the example
I provided. Can you check that the shapes of the two arrays are the same
in the code you execute? And if you can't find the origin of the error,
can you post a minimal example of code that reproduces the error?

	Also, as Anne mentioned (thanks for your answer, Anne!), the
weights keyword argument can be used with negative numbers, but also with
complex numbers (did you check it before posting?).
>>> a = np.arange(10)
>>> b = 1j*np.ones_like(a)
>>> np.histogram(a, a, weights=b)
(array([ 0.+1.j,  0.+1.j,  0.+1.j,  0.+1.j,  0.+1.j,  0.+1.j,  0.+1.j,
        0.+1.j,  0.+2.j]), array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9]))
So you don't need to separate real and imaginary parts, you can use
directly your complex array as weights.

	Cheers,

	Emmanuelle

> Cheers,

> Josh Lawrence
> Ph.D. Student
> Clemson University

> On Jun 12, 2009, at 2:28 PM, Emmanuelle Gouillart 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

> _______________________________________________
> SciPy-user mailing list
> SciPy-user at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user