On Fri, Feb 3, 2012 at 2:33 PM, <josef.p...@gmail.com> wrote: > On Fri, Feb 3, 2012 at 1:58 PM, santhu kumar <mesan...@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.27595761e-12, 0.00000000e+00, 3.63797881e-12], [ 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, <numpy-discussion-requ...@scipy.org> wrote: >>> >>> Send NumPy-Discussion mailing list submissions to >>> numpy-discussion@scipy.org >>> >>> To subscribe or unsubscribe via the World Wide Web, visit >>> http://mail.scipy.org/mailman/listinfo/numpy-discussion >>> or, via email, send a message with subject or body 'help' to >>> numpy-discussion-requ...@scipy.org >>> >>> You can reach the person managing the list at >>> numpy-discussion-ow...@scipy.org >>> >>> When replying, please edit your Subject line so it is more specific >>> than "Re: Contents of NumPy-Discussion digest..." >>> >>> >>> Today's Topics: >>> >>> 1. Re: Trick for fast (santhu kumar) >>> 2. Re: Trick for fast (josef.p...@gmail.com) >>> >>> >>> ---------------------------------------------------------------------- >>> >>> Message: 1 >>> Date: Fri, 3 Feb 2012 12:29:26 -0600 >>> From: santhu kumar <mesan...@gmail.com> >>> >>> Subject: Re: [Numpy-discussion] Trick for fast >>> To: numpy-discussion@scipy.org >>> Message-ID: >>> >>> <ca+7trsszzxplaoz_9pgu5e9mychtdqzmzrompzpfhrrkcst...@mail.gmail.com> >>> Content-Type: text/plain; charset="iso-8859-1" >>> >>> >>> On Fri, Feb 3, 2012 at 1:29 PM, santhu kumar <mesan...@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 >>> >>> >>> >>> >> >> _______________________________________________ >> NumPy-Discussion mailing list >> NumPy-Discussion@scipy.org >> http://mail.scipy.org/mailman/listinfo/numpy-discussion >> _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion