Add pyfaust demo scripts (ported from matlab: fft, runtimecmp).

wrapper/python/pyfaust/demo/fft/fig_speed_up_fourier.py

+
+from pylab import *
+
+
+def fig_speed_up_fourier():
+    from pyfaust.demo.fft.speed_up_fourier import nb_mults, log2_dims, dims, FFT_FAUST, \
+    FFT_FAUST_FULL, FFT_NATIVE, NUM_FFT_TYPES
+    import os.path, os, sys
+
+    output_dir = 'pyfaust_demo_figures'
+
+    input_dir = 'pyfaust_demo_output'
+    input_path = input_dir+os.sep+'speed_up_fourier.txt'
+    mult_times = loadtxt(input_path).reshape(nb_mults, len(dims), NUM_FFT_TYPES)
+    # NUM_FFT_TYPES == 3 is for FFT_FAUST, FFT_FAUST_FULL, FFT_NATIVE
+
+    mean_mult_times = squeeze(mult_times.mean(axis=0))
+
+    curve_thickness = 2
+    legend_loc = 'upper left'
+
+    plt.rcParams['figure.figsize'] = [12.0, 8]
+    subplot(211)
+
+    hold(True)
+    line_marker_types = [ 'ro-', 'bo-', 'go-', 'r+-', 'b+-', 'g+-' ]
+
+    title('Runtime Fourier A*x')
+    for i in range(0,NUM_FFT_TYPES):
+        semilogy(log2_dims, mean_mult_times[:,i],line_marker_types[i],
+                 lw=curve_thickness)
+        ylabel('Computed Time (sec)')
+    # legend(['dense', 'Faust', 'fft'], loc=legend_loc) # in first subplot
+    grid(True)
+    axes([log2_dims[0], log2_dims[-1], mean_mult_times.min(),
+          mean_mult_times.max()])
+    xticks(log2_dims)
+
+    subplot(212)
+    grid(True)
+    title('Speedup Fourier A*x')
+    for i in range(1,NUM_FFT_TYPES):
+        semilogy(log2_dims, mean_mult_times[:,0]/mean_mult_times[:,i],
+                 line_marker_types[i], lw=curve_thickness)
+        ylabel('Speedup Fourier A*x')
+        xlabel ('log(dim)')
+
+    figlegend(figure(1).get_axes()[0].get_lines(),['dense', 'faust',
+                                                   'fft'],loc='best')
+
+
+    if(not os.path.exists(output_dir)):
+        os.mkdir(output_dir)
+    savefig(output_dir+os.sep+'Fourier-RuntimeComp-Speedup.png',
+            dpi=200)
+    show()
+
+
+
+
+

wrapper/python/pyfaust/demo/fft/speed_up_fourier.py

+from __future__ import print_function
+from pylab import *
+import os,sys
+
+nb_mults = 500
+log2_dims = arange(6,12)
+dims = 2**log2_dims
+
+FFT_FAUST_FULL=0
+FFT_FAUST=1
+FFT_NATIVE=2
+
+NUM_FFT_TYPES=3
+
+def speed_up_fourier():
+    from pyfaust import Faust, FaustFactory
+    treshold = 10**-10
+
+    print('Speedup Fourier')
+    print("===============")
+    print('Generating data...')
+
+    fft_mats = []
+    fft_fausts = []
+
+    for k in range(0,len(dims)):
+        print('\rFFT dims processed: ', dims[0:k+1], end='')
+        F = FaustFactory.fourier(log2_dims[k])
+        fft_mat = F.toarray()#fft(eye(dims[k])) # or TODO: F.toarray() ?
+        fft_fausts += [ F ]
+        fft_mats += [ fft_mat ]
+    print()
+
+    print("Gathering computation times...")
+    from time import time, clock
+    if sys.platform == 'win32':
+        timer = clock
+    else:
+        timer = time
+
+    mult_times = ndarray(shape=(nb_mults, len(dims), 3))
+    # 3 for: FFT_FAUST_FULL, FFT_FAUST, FFT_NATIVE
+
+
+    for i in range(0,nb_mults):
+        print('\r#muls:',i+1,'/', nb_mults, end='')
+        for k in range(0,len(dims)):
+           dim = dims[k]
+           fft_mat = fft_mats[k]
+           fft_faust = fft_fausts[k]
+           ydense = empty(dim)
+           yfaust = empty(dim)
+           yfft = empty(dim)
+           ydense_trans = empty(dim)
+           yfaust_trans = empty(dim)
+           x = rand(dim,1)
+           t = timer()
+           ydense = fft_mat.dot(x)
+           mult_times[i,k,FFT_FAUST_FULL] = timer()-t
+           t = timer()
+           yfaust = fft_faust*x
+           mult_times[i,k,FFT_FAUST] = timer()-t
+           t = timer()
+           yfft = fft2(x)
+           mult_times[i,k,FFT_NATIVE] = timer()-t
+           if(norm(ydense-yfaust)>treshold):
+               raise Exception('Multiplication error: larger error than '
+                               'treshold for ydense.')
+
+           n1 = norm(fft(np.eye(dims[k]))-fft_faust.toarray())
+           assert(n1 < treshold)
+           if(norm(yfft-yfaust)/norm(yfft)>treshold):
+               print('\nerror:', norm(yfft-yfaust))
+               raise Exception('Multiplication error: larger error than '
+                               'treshold for yfaust')
+
+    print()
+
+
+    output_dir = 'pyfaust_demo_output'
+    import os.path
+    if(not os.path.exists(output_dir)):
+       os.mkdir(output_dir)
+    output_path = output_dir+os.sep+'speed_up_fourier.txt'
+    savetxt(output_path, mult_times.reshape(mult_times.size))
+
+    # test
+    # mult_times_r = loadtxt(output_path).reshape(mult_times.shape)
+    # assert(all(mult_times_r == mult_times))
+
+
+

wrapper/python/pyfaust/demo/runtimecmp/fig_runtime_comparison.py

+from pylab import *
+from pyfaust.demo.runtimecmp.runtime_comparison import rcgs, dims, nb_facts, nb_mults,\
+        dims_len, rcgs_len, nb_facts_len
+
+def fig_runtime_comparison():
+    import os.path, os
+    output_dir = 'pyfaust_demo_figures'
+    input_dir = 'pyfaust_demo_output'
+    input_file_existing = False
+
+    for matrix_or_vector in  ('vector', 'matrix'):
+        path_tfaust = input_dir+os.sep+'runtime_cmp_tfaust-'+matrix_or_vector+'.txt'
+        path_tdense = input_dir+os.sep+'runtime_cmp_tdense-'+matrix_or_vector+'.txt'
+        if(os.path.exists(path_tfaust) and os.path.exists(path_tdense)):
+            input_file_existing = True
+            break
+
+    if(not input_file_existing):
+        raise Exception("Input files don't exist, please run"
+                        " runtime_comparison.py first.")
+
+    if(not os.path.exists(output_dir)):
+        os.mkdir(output_dir)
+
+    tfaust = loadtxt(path_tfaust).reshape(nb_mults, dims_len, rcgs_len, nb_facts_len, 2)
+    tdense = loadtxt(path_tdense).reshape(nb_mults, dims_len, 2)
+
+
+    # average the mul. times according the number of runs
+    mean_tdense = squeeze(tdense.mean(axis=0))
+    mean_tfaust = squeeze(tfaust.mean(axis=0))
+
+    # avoid legend overlapping axis
+    plt.rcParams['figure.figsize'] = [12.0, 8]
+    #print("mean_tdense=", mean_tdense, mean_tdense.shape)
+
+# plot the time computed in logscale with a fixed number of factors in a given figure
+# in each figure we have in x-axis the log2(dimension) of the square matrix
+#                         in y-axis the time
+# all times for faust multiplication with different RCG (density)
+# and the times for dense matrix multiplication are plotted
+
+    curve_thickness = 2
+
+    # hold differents figures in the same box
+    ymin = min(mean_tdense.min(), mean_tfaust.min())
+    ymax = max(mean_tdense.max(), mean_tfaust.max())
+
+    legendy = [ 'Time (A*x)', 'Time (A.T*x)']
+
+    fig, ax = subplots(2,nb_facts_len, sharex=True, sharey=True)#, h*nb_facts_len+nf+1)
+    for h in arange(0,2):
+        for nf in range(0,nb_facts_len):
+
+            legend_curve = []
+            lines = []
+            for k in range(0, rcgs_len):
+                lines.append(*ax[h,nf].semilogy(log2(dims), mean_tfaust[:, k, nf, h],
+                         '-+', lw=curve_thickness))
+                legend_curve.append('Faust RCG '+str(rcgs[k]))
+                hold(True)
+
+            lines.append(*ax[h,nf].semilogy(log2(dims), squeeze(mean_tdense[:,h]), '-+', c=(0, .8, .8),
+                     lw=curve_thickness))
+
+            legend_curve.append('Dense ')
+
+            ax[h,nf].grid(True)
+            axes([log2(dims[0]), log2(dims[-1]), ymin, ymax])
+            if(h == 0):
+                ax[h,nf].set_title('#factors: '+str(nb_facts[nf]))
+            # legend for it axis (a bit heavy to read)
+            #ax[h,nf].legend(legend_curve)
+            #legend('best')
+            if(nf == 0):
+                if(h == 1):
+                    ax[h,nf].set_xlabel("log2(Dimension)")
+                ax[h,nf].set_ylabel(legendy[h])
+
+    figlegend(lines,legend_curve,loc='upper left')
+    # TODO: figure.Name matlab ?
+    constraint_str = 'sp'
+    # fig_name = 'Faust-'+matrix_or_vector+' multiplication '
+    #        '(constraint: 'constraint_str+')')
+    savefig(output_dir+os.sep+'RuntimeComp-'+matrix_or_vector+'_multiplication_constraint_'+constraint_str+'.png',
+           dpi=200)
+    #tight_layout()
+    #show()

wrapper/python/pyfaust/demo/runtimecmp/runtime_comparison.py

+from __future__ import print_function
+from pylab import *
+import os, sys
+
+rcgs = [2, 4, 8]
+dims = [128, 256, 512]
+nb_facts = [2, 4, 8]
+nb_mults = 500
+constraint = True # 'per_row' # per_col
+dims_len = len(dims)
+rcgs_len = len(rcgs)
+nb_facts_len = len(nb_facts)
+
+def runtime_comparison():
+    """
+    TODO
+    """
+    from pyfaust import Faust, FaustFactory
+
+    matrix_or_vector = 'vector'
+
+
+    fausts = ndarray(shape=(dims_len, rcgs_len, nb_facts_len), dtype=Faust)
+    dense_mats =  list()
+
+
+    # loading all the different Fausts and dense matrices
+    for j in range(0,dims_len):
+        dim = dims[j]
+        A = rand(dim, dim)
+        dense_mats.append(A)
+
+        for k in range(0,rcgs_len):
+            rcg = rcgs[k]
+
+            for l in range(0,nb_facts_len):
+                nf = nb_facts[l]
+                F = FaustFactory.rand(nf , dim, 1./(nf*rcg),
+                                      per_row=constraint, fac_type='sparse')
+                assert(F.rcg() == rcg)
+                fausts[j,k,l] = F
+                assert(F.shape == (dim,dim))
+
+    # time comparison
+    from time import time, clock
+    if sys.platform == 'win32':
+        timer = clock
+    else:
+        timer = time
+
+    tdense = ndarray(shape=(nb_mults, dims_len, 2))
+    tfaust = ndarray(shape=(nb_mults, dims_len, rcgs_len, nb_facts_len, 2))
+
+    for i in range(0, nb_mults):
+        print("\r\r #muls =",i+1,'/',nb_mults, end='')
+        for j in range(0,dims_len):
+            dim = dims[j]
+
+            if(matrix_or_vector == 'matrix'):
+                dim2 = dim # mul. by a square matrix
+            elif(matrix_or_vector == 'vector'):
+                dim2 = 1
+            else:
+                raise("matrix_or_vector string must be equal to matrix or"
+                      " vector")
+
+            for k in range(0,rcgs_len):
+                rcg = rcgs[k]
+
+                for l in range(0,nb_facts_len):
+                    nf = nb_facts[l]
+                    x = rand(dim,dim2)
+
+                    if(k == 0 and l == 0):
+                        A = dense_mats[j]
+                        t = timer()
+                        y = A.dot(x)
+                        tdense[i,j,0] = timer()-t
+                        t = timer()
+                        y_trans = A.T.dot(x)
+                        tdense[i,j,1] = timer()-t
+
+                    F = fausts[j,k,l]
+                    t = timer()
+                    yfaust = F*x
+                    tfaust[i,j,k,l,0] = timer()-t
+                    t = timer()
+                    yfaust_trans = F.T*x
+                    tfaust[i,j,k,l,1] = timer()-t
+
+    print()
+    output_dir = 'pyfaust_demo_output'
+    import os.path
+    if(not os.path.exists(output_dir)):
+       os.mkdir(output_dir)
+    path_tfaust = output_dir+os.sep+'runtime_cmp_tfaust-'+matrix_or_vector+'.txt'
+    path_tdense = output_dir+os.sep+'runtime_cmp_tdense-'+matrix_or_vector+'.txt'
+    savetxt(path_tdense, tdense.reshape(tdense.size))
+    savetxt(path_tfaust, tfaust.reshape(tfaust.size))
+
+    # test
+    #tfaust_r = loadtxt(path_tfaust).reshape(tfaust.shape)
+    #assert(all(tfaust_r == tfaust))
+    #tdense_r = loadtxt(path_tdense).reshape(tdense.shape)
+    #assert(all(tdense_r == tdense))
+
+
+