Hi Josef, I am unclear on what you want to say, but all I am doing in the code is getting inertia tensor for a bunch of particle masses. (http://en.wikipedia.org/wiki/Moment_of_inertia#Moment_of_inertia_tensor) So the diagonals are not actually zeros but would have z^2 + y^2 .. The reason which I said 3secs could be misunderstood .. This code is called many times over a loop and the bulk time is taken in computing this inertial tensor. After code change, the entire loop finishes off in 3 ses. Thanks for alertness, Santhosh On Fri, Feb 3, 2012 at 12:47 PM, <numpydiscussionrequest@scipy.org> wrote:
Send NumPyDiscussion mailing list submissions to numpydiscussion@scipy.org
To subscribe or unsubscribe via the World Wide Web, visit http://mail.scipy.org/mailman/listinfo/numpydiscussion or, via email, send a message with subject or body 'help' to numpydiscussionrequest@scipy.org
You can reach the person managing the list at numpydiscussionowner@scipy.org
When replying, please edit your Subject line so it is more specific than "Re: Contents of NumPyDiscussion digest..."
Today's Topics:
1. Re: Trick for fast (santhu kumar) 2. Re: Trick for fast (josef.pktd@gmail.com)

Message: 1 Date: Fri, 3 Feb 2012 12:29:26 0600 From: santhu kumar <mesanthu@gmail.com> Subject: Re: [Numpydiscussion] Trick for fast To: numpydiscussion@scipy.org MessageID: <CA+7TRsszzxPLAOz_9pgU5e9mYcHTDQZMzRoMPzpfhrRkcStZEQ@mail.gmail.com
ContentType: text/plain; charset="iso88591"
On Fri, Feb 3, 2012 at 1:29 PM, santhu kumar <mesanthu@gmail.com> wrote:
Hello all,
Thanks for lovely solutions. I have sat on it for some time and wrote it myself :
n =x.shape[0] ea = np.array([1,0,0,0,1,0,0,0,1]) inert = ((np.tile(ea,(n,1))*((x*x).sum(axis=1)[:,np.newaxis]) 
np.hstack([x*x[:,0][:,np.newaxis],x*x[:,1][:,np.newaxis],x*x[:,2][:,np.newaxis]]))*mass[:,np.newaxis]).sum(axis=0)
inert.shape = 3,3
Does the trick and reduces the time from over 45 secs to 3 secs. I do want to try einsum but my numpy is little old and it does not have it.
Thanks Sebastian (it was tricky to understand your code for me) and Josef (clean).
Isn't the entire substraction of the first term just to set the diagonal of the result to zero.
It looks to me now just like the weighted dot product and setting the diagonal to zero. That shouldn't take 3 secs unless you actual dimensions are huge.
Josef
On Fri, Feb 3, 2012 at 1:58 PM, santhu kumar <mesanthu@gmail.com> wrote:
Hi Josef,
I am unclear on what you want to say, but all I am doing in the code is getting inertia tensor for a bunch of particle masses. (http://en.wikipedia.org/wiki/Moment_of_inertia#Moment_of_inertia_tensor)
So the diagonals are not actually zeros but would have z^2 + y^2 .. The reason which I said 3secs could be misunderstood .. This code is called many times over a loop and the bulk time is taken in computing this inertial tensor. After code change, the entire loop finishes off in 3 ses.
ok, I still had a python sum instead of the numpy sum in there (I really really like namespaces :)
np.sum(ri*ri) 5549.0 sum(ri*ri) array([ 1764., 1849., 1936.])
this might match better. print (mass * x**2).sum() * np.eye(3)  np.dot(x.T, mass*x) or res =  np.dot(x.T, mass*x) res[np.arange(3), np.arange(3)] += (mass * x**2).sum() print res Josef
Thanks for alertness, Santhosh
On Fri, Feb 3, 2012 at 12:47 PM, <numpydiscussionrequest@scipy.org> wrote:
Send NumPyDiscussion mailing list submissions to numpydiscussion@scipy.org
To subscribe or unsubscribe via the World Wide Web, visit http://mail.scipy.org/mailman/listinfo/numpydiscussion or, via email, send a message with subject or body 'help' to numpydiscussionrequest@scipy.org
You can reach the person managing the list at numpydiscussionowner@scipy.org
When replying, please edit your Subject line so it is more specific than "Re: Contents of NumPyDiscussion digest..."
Today's Topics:
1. Re: Trick for fast (santhu kumar) 2. Re: Trick for fast (josef.pktd@gmail.com)

Message: 1 Date: Fri, 3 Feb 2012 12:29:26 0600 From: santhu kumar <mesanthu@gmail.com>
Subject: Re: [Numpydiscussion] Trick for fast To: numpydiscussion@scipy.org MessageID:
<CA+7TRsszzxPLAOz_9pgU5e9mYcHTDQZMzRoMPzpfhrRkcStZEQ@mail.gmail.com> ContentType: text/plain; charset="iso88591"
On Fri, Feb 3, 2012 at 1:29 PM, santhu kumar <mesanthu@gmail.com> wrote:
Hello all,
Thanks for lovely solutions. I have sat on it for some time and wrote it myself :
n =x.shape[0] ea = np.array([1,0,0,0,1,0,0,0,1]) inert = ((np.tile(ea,(n,1))*((x*x).sum(axis=1)[:,np.newaxis]) 
np.hstack([x*x[:,0][:,np.newaxis],x*x[:,1][:,np.newaxis],x*x[:,2][:,np.newaxis]]))*mass[:,np.newaxis]).sum(axis=0) inert.shape = 3,3
Does the trick and reduces the time from over 45 secs to 3 secs. I do want to try einsum but my numpy is little old and it does not have it.
Thanks Sebastian (it was tricky to understand your code for me) and Josef (clean).
Isn't the entire substraction of the first term just to set the diagonal of the result to zero.
It looks to me now just like the weighted dot product and setting the diagonal to zero. That shouldn't take 3 secs unless you actual dimensions are huge.
Josef
_______________________________________________ NumPyDiscussion mailing list NumPyDiscussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpydiscussion
On Fri, Feb 3, 2012 at 2:33 PM, <josef.pktd@gmail.com> wrote:
On Fri, Feb 3, 2012 at 1:58 PM, santhu kumar <mesanthu@gmail.com> wrote:
Hi Josef,
I am unclear on what you want to say, but all I am doing in the code is getting inertia tensor for a bunch of particle masses. (http://en.wikipedia.org/wiki/Moment_of_inertia#Moment_of_inertia_tensor)
So the diagonals are not actually zeros but would have z^2 + y^2 .. The reason which I said 3secs could be misunderstood .. This code is called many times over a loop and the bulk time is taken in computing this inertial tensor. After code change, the entire loop finishes off in 3 ses.
ok, I still had a python sum instead of the numpy sum in there (I really really like namespaces :)
np.sum(ri*ri) 5549.0 sum(ri*ri) array([ 1764., 1849., 1936.])
this might match better.
print (mass * x**2).sum() * np.eye(3)  np.dot(x.T, mass*x)
or res =  np.dot(x.T, mass*x) res[np.arange(3), np.arange(3)] += (mass * x**2).sum() print res
if my interpretation is now correct, and because I have seen many traces lately
res =  np.dot(x.T, mass*x) res[np.arange(3), np.arange(3)] = np.trace(res) res array([[ 35752.5, 16755. , 17287.5], [16755. , 34665. , 17865. ], [17287.5, 17865. , 33532.5]])
res  inert array([[ 7.27595761e12, 0.00000000e+00, 3.63797881e12], [ 0.00000000e+00, 0.00000000e+00, 0.00000000e+00], [ 0.00000000e+00, 0.00000000e+00, 0.00000000e+00]])
Josef
Josef
Thanks for alertness, Santhosh
On Fri, Feb 3, 2012 at 12:47 PM, <numpydiscussionrequest@scipy.org> wrote:
Send NumPyDiscussion mailing list submissions to numpydiscussion@scipy.org
To subscribe or unsubscribe via the World Wide Web, visit http://mail.scipy.org/mailman/listinfo/numpydiscussion or, via email, send a message with subject or body 'help' to numpydiscussionrequest@scipy.org
You can reach the person managing the list at numpydiscussionowner@scipy.org
When replying, please edit your Subject line so it is more specific than "Re: Contents of NumPyDiscussion digest..."
Today's Topics:
1. Re: Trick for fast (santhu kumar) 2. Re: Trick for fast (josef.pktd@gmail.com)

Message: 1 Date: Fri, 3 Feb 2012 12:29:26 0600 From: santhu kumar <mesanthu@gmail.com>
Subject: Re: [Numpydiscussion] Trick for fast To: numpydiscussion@scipy.org MessageID:
<CA+7TRsszzxPLAOz_9pgU5e9mYcHTDQZMzRoMPzpfhrRkcStZEQ@mail.gmail.com> ContentType: text/plain; charset="iso88591"
On Fri, Feb 3, 2012 at 1:29 PM, santhu kumar <mesanthu@gmail.com> wrote:
Hello all,
Thanks for lovely solutions. I have sat on it for some time and wrote it myself :
n =x.shape[0] ea = np.array([1,0,0,0,1,0,0,0,1]) inert = ((np.tile(ea,(n,1))*((x*x).sum(axis=1)[:,np.newaxis]) 
np.hstack([x*x[:,0][:,np.newaxis],x*x[:,1][:,np.newaxis],x*x[:,2][:,np.newaxis]]))*mass[:,np.newaxis]).sum(axis=0) inert.shape = 3,3
Does the trick and reduces the time from over 45 secs to 3 secs. I do want to try einsum but my numpy is little old and it does not have it.
Thanks Sebastian (it was tricky to understand your code for me) and Josef (clean).
Isn't the entire substraction of the first term just to set the diagonal of the result to zero.
It looks to me now just like the weighted dot product and setting the diagonal to zero. That shouldn't take 3 secs unless you actual dimensions are huge.
Josef
_______________________________________________ NumPyDiscussion mailing list NumPyDiscussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpydiscussion
On Fri, Feb 3, 2012 at 4:49 PM, Alan G Isaac <alan.isaac@gmail.com> wrote:
On 2/3/2012 3:37 PM, josef.pktd@gmail.com wrote:
res =  np.dot(x.T, mass*x) res[np.arange(3), np.arange(3)] = np.trace(res)
Nice! Get some speed gain with slicing:
res =  np.dot(x.T, mass*x) res.flat[slice(0,None,4)] = np.trace(res)
Actually, I thought about the most readable version using diag_indices
di = np.diag_indices(3) res =  np.dot(x.T, mass*x) res[di] = np.trace(res)
(but I thought I said already enough) Josef
Alan
_______________________________________________ NumPyDiscussion mailing list NumPyDiscussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpydiscussion
participants (3)

Alan G Isaac

josef.pktd＠gmail.com

santhu kumar