Re: r23935 - /branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/ -- June 13, 2014 - 14:44

Hi,

I've been doing this for a long time.  In the last 2 years I have made
the calculation for the most complex frame order model - the theory
I'm developing and will hopefully publish soon - drop from 1,000,000
years down to 1 day!  I also spend a lot of time during my PhD
optimising the model-free target functions, gradients, and Hessians
for the fastest speed possible.  The PDB reading is another area that
I've massively speed up, though speeding up the reading text files is
a completely different game.  I've never added my testing scripts
before and thought, after you added yours, that it is a great
demonstration.  This is also very useful for the future as numpy is
optimised more.  Some operations that are slow today might be faster
than all others with some future numpy version.  I've been observing
these speed changes since the Numeric to numpy transition that I had
to guide relax through, many, many years ago.

The timeit module is different to profiling, as it runs the code X
times rather than once, so the fluctuations in the timings are better
averaged compared to the profiling.  But the profiling does not suffer
from the CPU being used by another process whereas timeit does.  So it
really depends on the system and how long one iteration takes.  I like
the timeit module in this case as everything is run one after the
other within one script, so that CPU usage fluctuations are not a
significant effect.  And the summary is much better with a single
timing in seconds printed out.  So both profiling and timeit are
useful.

Regards,

Edward





On 13 June 2014 14:18, Troels Emtekær Linnet <tlinnet@xxxxxxxxxxxxx> wrote:
> Seems like I am not the only one having fun with higher dimensions. ;-)
>
> is timeit better than profiling?
>
> How about the model_free analysis?
>
> Is that in for loops?
>
> best
> Troels
>
>
> 2014-06-13 14:03 GMT+02:00 <edward@xxxxxxxxxxxxx>:
>
> > Author: bugman
> > Date: Fri Jun 13 14:03:44 2014
> > New Revision: 23935
> >
> > URL: http://svn.gna.org/viewcvs/relax?rev=23935&view=rev
> > Log:
> > Added a script for timing different ways to calculate PCSs and RDCs for
> > multiple vectors.
> >
> > This uses the timeit module rather than profile to demonstrate the speed
> > of 7 different ways to
> > calculate the RDCs or PCSs for an array of vectors using numpy.  In the
> > frame order analysis, this
> > is the bottleneck for the quasi-random numerical integration of the PCS.
> >
> > The log file shows a potential 1 order of magnitude speed up between the
> > 1st technique, which is
> > currently used in the frame order analysis, and the 7th and last
> > technique.  The first technique
> > loops over each vector, calculating the PCS.  The last expands the PCS/RDC
> > equation of the
> > projection of the vector into the alignment tensor, and calculates all
> > PCSs simultaneously.
> >
> >
> > Added:
> >
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/
> >
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.log
> >
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.py
> >
> > Added:
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.log
> > URL:
> > http://svn.gna.org/viewcvs/relax/branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.log?rev=23935&view=auto
> >
> > ==============================================================================
> > ---
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.log
> >      (added)
> > +++
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.log
> >      Fri Jun 13 14:03:44 2014
> > @@ -0,0 +1,38 @@
> > +Original vectors:
> > +[[ 1.  2.  3.]
> > + [ 2.  2.  2.]]
> > +Shape: (200, 3)
> > +
> > +Tensor:
> > +[[ 1.42219822 -7.07796212 -6.01619494]
> > + [-7.07796212 -1.44543002  2.02008007]
> > + [-6.01619494  2.02008007  0.02323179]]
> > +
> > +
> > +1st projection - element by element r[i].A.r[i].
> > +Proj: [-44.31849296 -88.59261589]
> > +Timing (s): 3.90826106071
> > +
> > +2nd projection - diag of double tensordot.
> > +Proj: [-44.31849296 -88.59261589]
> > +Timing (s): 2.68392705917
> > +
> > +3rd projection - diag of double tensordot, no transpose.
> > +Proj: [-44.31849296 -88.59261589]
> > +Timing (s): 2.58908486366
> > +
> > +4th projection - mixed tensordot() and per-vector dot().
> > +Proj: [-44.31849296 -88.59261589]
> > +Timing (s): 5.47785592079
> > +
> > +5th projection - expansion and sum.
> > +Proj: [-44.31849296 -88.59261589]
> > +Timing (s): 10.7986190319
> > +
> > +6th projection - expansion.
> > +Proj: [-44.31849296 -88.59261589]
> > +Timing (s): 0.54491686821
> > +
> > +7th projection - expansion with pre-calculation.
> > +Proj: [-44.31849296 -88.59261589]
> > +Timing (s): 0.395635128021
> >
> > Added:
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.py
> > URL:
> > http://svn.gna.org/viewcvs/relax/branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.py?rev=23935&view=auto
> >
> > ==============================================================================
> > ---
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.py
> >       (added)
> > +++
> > branches/frame_order_cleanup/test_suite/shared_data/frame_order/timings/tensor_projections.py
> >       Fri Jun 13 14:03:44 2014
> > @@ -0,0 +1,127 @@
> > +# Python module imports.
> > +from numpy import *
> > +from numpy.linalg import norm
> > +from os import pardir, sep
> > +import sys
> > +from time import sleep
> > +from timeit import timeit
> > +
> > +# Modify the system path.
> > +sys.path.append(pardir+sep+pardir+sep+pardir+sep+pardir+sep)
> > +
> > +# relax module imports.
> > +from lib.geometry.rotations import euler_to_R_zyz
> > +
> > +
> > +def proj1(vect, A, N=1, verb=True):
> > +    d = zeros(len(vect), float64)
> > +    for i in range(N):
> > +        for j in xrange(len(vect)):
> > +            d[j] = dot(vect[j], dot(A, vect[j]))
> > +    if verb:
> > +        print("\n1st projection - element by element r[i].A.r[i].")
> > +        print("Proj: %s" % d[:2])
> > +
> > +
> > +def proj2(vect, A, N=1, verb=True):
> > +    for i in range(N):
> > +        d = diagonal(tensordot(vect, tensordot(A, transpose(vect),
> > axes=1), axes=1))
> > +    if verb:
> > +        print("\n2nd projection - diag of double tensordot.")
> > +        print("Proj: %s" % d[:2])
> > +
> > +
> > +def proj3(vect, A, N=1, verb=True):
> > +    for i in range(N):
> > +        d = diagonal(tensordot(tensordot(A, vect, axes=([0],[1])), vect,
> > axes=([0],[1])))
> > +    if verb:
> > +        print("\n3rd projection - diag of double tensordot, no
> > transpose.")
> > +        print("Proj: %s" % d[:2])
> > +
> > +
> > +def proj4(vect, A, N=1, verb=True):
> > +    d = zeros(len(vect), float64)
> > +    for i in range(N):
> > +        a = tensordot(A, vect, axes=([0],[1]))
> > +        for j in range(len(vect)):
> > +            d[j] = dot(vect[j], a[:,j])
> > +    if verb:
> > +        print("\n4th projection - mixed tensordot() and per-vector
> > dot().")
> > +        print("Proj: %s" % d[:2])
> > +
> > +
> > +def proj5(vect, A, N=1, verb=True):
> > +    d = zeros(len(vect), float64)
> > +    for i in range(N):
> > +        vect2 = vect**2
> > +        double_vect = 2.0 * vect
> > +        for j in xrange(len(vect)):
> > +            d[j] = vect2[j][0]*A[0, 0] + vect2[j][1]*A[1, 1] +
> > vect2[j][2]*(A[2, 2]) + double_vect[j][0]*vect[j][1]*A[0, 1] +
> > double_vect[j][0]*vect[j][2]*A[0, 2] + double_vect[j][1]*vect[j][2]*A[1, 2]
> > +    if verb:
> > +        print("\n5th projection - expansion and sum.")
> > +        print("Proj: %s" % d[:2])
> > +
> > +
> > +def proj6(vect, A, N=1, verb=True):
> > +    d = zeros(len(vect), float64)
> > +    for i in range(N):
> > +        d = vect[:,0]**2*A[0, 0] + vect[:,1]**2*A[1, 1] +
> > vect[:,2]**2*(A[2, 2]) + 2.0*vect[:,0]*vect[:,1]*A[0, 1] +
> > 2.0*vect[:,0]*vect[:,2]*A[0, 2] + 2.0*vect[:,1]*vect[:,2]*A[1, 2]
> > +    if verb:
> > +        print("\n6th projection - expansion.")
> > +        print("Proj: %s" % d[:2])
> > +
> > +
> > +def proj7(vect, A, N=1, verb=True):
> > +    d = zeros(len(vect), float64)
> > +    for i in range(N):
> > +        vect2 = vect**2
> > +        double_vect = 2.0 * vect
> > +        d = vect2[:,0]*A[0, 0] + vect2[:,1]*A[1, 1] + vect2[:,2]*(A[2,
> > 2]) + double_vect[:,0]*vect[:,1]*A[0, 1] + double_vect[:,0]*vect[:,2]*A[0,
> > 2] + double_vect[:,1]*vect[:,2]*A[1, 2]
> > +    if verb:
> > +        print("\n7th projection - expansion with pre-calculation.")
> > +        print("Proj: %s" % d[:2])
> > +
> > +
> > +# Some 200 vectors.
> > +vect = array([[1, 2, 3], [2, 2, 2]], float64)
> > +vect = tile(vect, (100, 1))
> > +if __name__ == '__main__':
> > +    print("Original vectors:\n%s\nShape: %s" % (vect[:2], vect.shape))
> > +
> > +# Init the alignment tensor.
> > +A = zeros((3, 3), float64)
> > +A_5D = [1.42219822168827662867e-04, -1.44543001566521341940e-04,
> > -7.07796211648713973798e-04, -6.01619494082773244303e-04,
> > 2.02008007072950861996e-04]
> > +A[0, 0] = A_5D[0]
> > +A[1, 1] = A_5D[1]
> > +A[2, 2] = -A_5D[0] -A_5D[1]
> > +A[0, 1] = A[1, 0] = A_5D[2]
> > +A[0, 2] = A[2, 0] = A_5D[3]
> > +A[1, 2] = A[2, 1] = A_5D[4]
> > +A *= 1e4
> > +if __name__ == '__main__':
> > +    print("\nTensor:\n%s\n" % A)
> > +
> > +# Projections.
> > +N = 100
> > +M = 100
> > +if __name__ == '__main__':
> > +    proj1(vect=vect, A=A, N=1, verb=True)
> > +    print("Timing (s): %s" % timeit("proj1(vect=vect, A=A, N=N,
> > verb=False)", setup="from tensor_projections import proj1, vect, A, N",
> > number=M))
> > +
> > +    proj2(vect=vect, A=A, N=1, verb=True)
> > +    print("Timing (s): %s" % timeit("proj2(vect=vect, A=A, N=N,
> > verb=False)", setup="from tensor_projections import proj2, vect, A, N",
> > number=M))
> > +
> > +    proj3(vect=vect, A=A, N=1, verb=True)
> > +    print("Timing (s): %s" % timeit("proj3(vect=vect, A=A, N=N,
> > verb=False)", setup="from tensor_projections import proj3, vect, A, N",
> > number=M))
> > +
> > +    proj4(vect=vect, A=A, N=1, verb=True)
> > +    print("Timing (s): %s" % timeit("proj4(vect=vect, A=A, N=N,
> > verb=False)", setup="from tensor_projections import proj4, vect, A, N",
> > number=M))
> > +
> > +    proj5(vect=vect, A=A, N=1, verb=True)
> > +    print("Timing (s): %s" % timeit("proj5(vect=vect, A=A, N=N,
> > verb=False)", setup="from tensor_projections import proj5, vect, A, N",
> > number=M))
> > +
> > +    proj6(vect=vect, A=A, N=1, verb=True)
> > +    print("Timing (s): %s" % timeit("proj6(vect=vect, A=A, N=N,
> > verb=False)", setup="from tensor_projections import proj6, vect, A, N",
> > number=M))
> > +
> > +    proj7(vect=vect, A=A, N=1, verb=True)
> > +    print("Timing (s): %s" % timeit("proj7(vect=vect, A=A, N=N,
> > verb=False)", setup="from tensor_projections import proj7, vect, A, N",
> > number=M))
> >
> >
> > _______________________________________________
> > relax (http://www.nmr-relax.com)
> >
> > This is the relax-commits mailing list
> > relax-commits@xxxxxxx
> >
> > To unsubscribe from this list, get a password
> > reminder, or change your subscription options,
> > visit the list information page at
> > https://mail.gna.org/listinfo/relax-commits
> _______________________________________________
> relax (http://www.nmr-relax.com)
>
> This is the relax-devel mailing list
> relax-devel@xxxxxxx
>
> To unsubscribe from this list, get a password
> reminder, or change your subscription options,
> visit the list information page at
> https://mail.gna.org/listinfo/relax-devel