def main():

    N = 1e8 // 2
    print("Data size", N)

    targets = ["cpu", "parallel"]

    # run just one target if is specified in the argument
    for t in targets:
        if t in sys.argv[1:]:
            targets = [t]
            break

    for target in targets:
        print("== Target", target)
        vect_discriminant = vectorize(["f4(f4, f4, f4)", "f8(f8, f8, f8)"], target=target)(discriminant)

        A, B, C = generate_input(N, dtype=np.float32)
        D = np.empty(A.shape, dtype=A.dtype)

        ts = time()
        D = vect_discriminant(A, B, C)
        te = time()

        total_time = te - ts

        print("Execution time %.4f" % total_time)
        print("Throughput %.4f" % (N / total_time))

        if "-verify" in sys.argv[1:]:
            check_answer(D, A, B, C)

def run_div_test(fld, exact, title='', show=False, ignore_inexact=False):
    t0 = time()
    result_numexpr = viscid.div(fld, preferred="numexpr", only=False)
    t1 = time()
    logger.info("numexpr magnitude runtime: %g", t1 - t0)

    result_diff = viscid.diff(result_numexpr, exact)['x=1:-1, y=1:-1, z=1:-1']
    if not ignore_inexact and not (result_diff.data < 5e-5).all():
        logger.warning("numexpr result is far from the exact result")
    logger.info("min/max(abs(numexpr - exact)): %g / %g",
                np.min(result_diff.data), np.max(result_diff.data))

    planes = ["y=0j", "z=0j"]
    nrows = 2
    ncols = len(planes)
    _, axes = plt.subplots(nrows, ncols, squeeze=False)

    for i, p in enumerate(planes):
        vlt.plot(result_numexpr, p, ax=axes[0, i], show=False)
        vlt.plot(result_diff, p, ax=axes[1, i], show=False)

    plt.suptitle(title)
    vlt.auto_adjust_subplots(subplot_params=dict(top=0.9))

    plt.savefig(next_plot_fname(__file__))
    if show:
        vlt.mplshow()

def timeit(f, *args, **kwargs):
    t0 = time()
    ret = f(*args, **kwargs)
    t1 = time()

    print("Took {0:.03g} secs.".format(t1 - t0))
    return ret

    def test_func(self):
        A = np.array(np.random.random((n, n)), dtype=np.float32)
        B = np.array(np.random.random((n, n)), dtype=np.float32)
        C = np.empty_like(A)

        print("N = %d x %d" % (n, n))

        s = time()
        stream = cuda.stream()
        with stream.auto_synchronize():
            dA = cuda.to_device(A, stream)
            dB = cuda.to_device(B, stream)
            dC = cuda.to_device(C, stream)
            cu_square_matrix_mul[(bpg, bpg), (tpb, tpb), stream](dA, dB, dC)
            dC.copy_to_host(C, stream)

        e = time()
        tcuda = e - s

        # Host compute
        Amat = np.matrix(A)
        Bmat = np.matrix(B)

        s = time()
        Cans = Amat * Bmat
        e = time()
        tcpu = e - s

        print('cpu:  %f' % tcpu)
        print('cuda: %f' % tcuda)
        print('cuda speedup: %.2fx' % (tcpu / tcuda))

        # Check result
        self.assertTrue(np.allclose(C, Cans))

 def func2(num):
     s = time()
     num = sp.sympify(num)
     res = num.is_Symbol
     e = time()
     print e-s
     return res

    def _bin_data(self, X, rng, is_training_data):
        """Bin data X.

        If is_training_data, then set the bin_mapper_ attribute.
        Else, the binned data is converted to a C-contiguous array.
        """

        description = 'training' if is_training_data else 'validation'
        if self.verbose:
            print("Binning {:.3f} GB of {} data: ".format(
                X.nbytes / 1e9, description), end="", flush=True)
        tic = time()
        if is_training_data:
            X_binned = self.bin_mapper_.fit_transform(X)  # F-aligned array
        else:
            X_binned = self.bin_mapper_.transform(X)  # F-aligned array
            # We convert the array to C-contiguous since predicting is faster
            # with this layout (training is faster on F-arrays though)
            X_binned = np.ascontiguousarray(X_binned)
        toc = time()
        if self.verbose:
            duration = toc - tic
            print("{:.3f} s".format(duration))

        return X_binned

def main():
    targets = ['cpu', 'parallel']

    # run just one target if is specified in the argument
    for t in targets:
        if t in sys.argv[1:]:
            targets = [t]
            break

    for target in targets:
        print('== Target', target)
        vect_sum = vectorize(['f4(f4, f4)', 'f8(f8, f8)'],
                             target=target)(sum)

        A = np.random.random(N).astype(np.float32)
        B = np.random.random(N).astype(np.float32)
        assert A.shape == B.shape
        assert A.dtype ==  B.dtype
        assert len(A.shape) == 1

        D = np.empty(A.shape, dtype=A.dtype)

        print('Data size', N)

        ts = time()
        D = vect_sum(A, B)
        te = time()

        total_time = (te - ts)

        print('Execution time %.4f' % total_time)
        print('Throughput %.4f' % (N / total_time))

        if '-verify' in sys.argv[1:]:
            check_answer(D, A, B, C)

        def test(ty):
            print("Test %s" % ty)
            data = np.array(np.random.random(1e6 + 1), dtype=ty)

            ts = time()
            stream = cuda.stream()
            device_data = cuda.to_device(data, stream)
            dresult = cuda_ufunc(device_data, device_data, stream=stream)
            result = dresult.copy_to_host()
            stream.synchronize()
            tnumba = time() - ts

            ts = time()
            gold = np_ufunc(data, data)
            tnumpy = time() - ts

            print("Numpy time: %fs" % tnumpy)
            print("Numba time: %fs" % tnumba)

            if tnumba < tnumpy:
                print("Numba is FASTER by %fx" % (tnumpy / tnumba))
            else:
                print("Numba is SLOWER by %fx" % (tnumba / tnumpy))

            self.assertTrue(np.allclose(gold, result), (gold, result))

def timeit(message, display=True):
    """Context to time an execution."""
    start = time()
    yield
    if not display:
        return
    print("{}: {:.3f} s".format(message, time() - start))

 def run(self):
     print('Running part 5')
 
     filename = './' + self.out_dir + '/time.txt'
     with open(filename, 'w') as text_file:
         
         t0 = time()
         self.nn_pca_cluster_wine()
         text_file.write('nn_pca_wine: %0.3f seconds\n' % (time() - t0))
         
         t0 = time()
         self.nn_ica_cluster_wine()
         text_file.write('nn_ica_wine: %0.3f seconds\n' % (time() - t0))
         
         t0 = time()
         self.nn_rp_cluster_wine()
         text_file.write('nn_rp_wine: %0.3f seconds\n' % (time() - t0))
         
         t0 = time()
         self.nn_lda_cluster_wine()
         text_file.write('nn_lda_wine: %0.3f seconds\n' % (time() - t0))
         
         t0 = time()
         self.nn_wine_orig()
         text_file.write('nn_wine_orig: %0.3f seconds\n' % (time() - t0))

    def test_gufunc(self):

        @guvectorize([void(float32[:, :], float32[:, :], float32[:, :])],
                     '(m,n),(n,p)->(m,p)',
                     target='cuda')
        def matmulcore(A, B, C):
            m, n = A.shape
            n, p = B.shape
            for i in range(m):
                for j in range(p):
                    C[i, j] = 0
                    for k in range(n):
                        C[i, j] += A[i, k] * B[k, j]

        gufunc = matmulcore
        gufunc.max_blocksize = 512

        matrix_ct = 1001 # an odd number to test thread/block division in CUDA
        A = np.arange(matrix_ct * 2 * 4, dtype=np.float32).reshape(matrix_ct, 2,
                                                                   4)
        B = np.arange(matrix_ct * 4 * 5, dtype=np.float32).reshape(matrix_ct, 4,
                                                                   5)

        ts = time()
        C = gufunc(A, B)
        tcuda = time() - ts

        ts = time()
        Gold = ut.matrix_multiply(A, B)
        tcpu = time() - ts

        non_stream_speedups.append(tcpu / tcuda)

        self.assertTrue(np.allclose(C, Gold))

 def wrapper(self, *args, **kwargs):
     # Open link
     if self.link is None:
         self.open()
     # Init time
     start = time()
     no_control = self.callback_timeout >= self.instrument_timeout
     # Loop over timeouts
     while True:
         try:
             # Run the function
             result = func(self, *args, **kwargs)
         except Vxi11Exception as exc:
             # Time control
             no_timeout = exc.err != ERR_IO_TIMEOUT
             expired = time() > start + self.instrument_timeout
             # Reraise exception
             if no_control or no_timeout or expired:
                 raise
             # Callback with exc
             if self.callback:
                 self.callback(exc)
         else:
             # Callback without exc
             if self.callback:
                 self.callback(None)
             # Return
             return result

def main():

    N = 1e+8 // 2
    print('Data size', N)

    targets = ['cpu', 'parallel']
    
    # run just one target if is specified in the argument
    for t in targets:
        if t in sys.argv[1:]:
            targets = [t]
            break

    for target in targets:
        print('== Target', target)
        vect_discriminant = vectorize([f4(f4, f4, f4), f8(f8, f8, f8)],
                                    target=target)(discriminant)

        A, B, C = generate_input(N, dtype=np.float32)
        D = np.empty(A.shape, dtype=A.dtype)

        ts = time()
        D = vect_discriminant(A, B, C)
        te = time()

        total_time = (te - ts)

        print('Execution time %.4f' % total_time)
        print('Throughput %.4f' % (N / total_time))



        if '-verify' in sys.argv[1:]:
            check_answer(D, A, B, C)

    def test_gufunc_stream(self):
        #cuda.driver.flush_pending_free()
        matrix_ct = 1001 # an odd number to test thread/block division in CUDA
        A = np.arange(matrix_ct * 2 * 4, dtype=np.float32).reshape(matrix_ct, 2,
                                                                   4)
        B = np.arange(matrix_ct * 4 * 5, dtype=np.float32).reshape(matrix_ct, 4,
                                                                   5)

        ts = time()
        stream = cuda.stream()
        dA = cuda.to_device(A, stream)
        dB = cuda.to_device(B, stream)

        dC = cuda.device_array(shape=(1001, 2, 5), dtype=A.dtype, stream=stream)
        dC = gufunc(dA, dB, out=dC, stream=stream)
        C = dC.copy_to_host(stream=stream)
        stream.synchronize()

        tcuda = time() - ts

        ts = time()
        Gold = ut.matrix_multiply(A, B)
        tcpu = time() - ts

        stream_speedups.append(tcpu / tcuda)

        self.assertTrue(np.allclose(C, Gold))

    def test_func(self):
        np.random.seed(42)
        A = np.array(np.random.random((n, n)), dtype=np.float32)
        B = np.array(np.random.random((n, n)), dtype=np.float32)
        C = np.empty_like(A)

        s = time()
        stream = cuda.stream()
        with stream.auto_synchronize():
            dA = cuda.to_device(A, stream)
            dB = cuda.to_device(B, stream)
            dC = cuda.to_device(C, stream)
            cu_square_matrix_mul[(bpg, bpg), (tpb, tpb), stream](dA, dB, dC)
            dC.copy_to_host(C, stream)

        e = time()
        tcuda = e - s

        # Host compute
        s = time()
        Cans = np.dot(A, B)
        e = time()
        tcpu = e - s

        # Check result
        np.testing.assert_allclose(C, Cans, rtol=1e-5)

def run_mag_test(fld, title="", show=False):
    vx, vy, vz = fld.component_views()  # pylint: disable=W0612
    vx, vy, vz = fld.component_fields()

    try:
        t0 = time()
        mag_ne = viscid.magnitude(fld, preferred="numexpr", only=False)
        t1 = time()
        logger.info("numexpr mag runtime: %g", t1 - t0)
    except viscid.verror.BackendNotFound:
        xfail("Numexpr is not installed")

    planes = ["z=0", "y=0"]
    nrows = 4
    ncols = len(planes)

    _, axes = plt.subplots(nrows, ncols, sharex=True, sharey=True, squeeze=False)

    for ind, p in enumerate(planes):
        vlt.plot(vx, p, ax=axes[0, ind], show=False)
        vlt.plot(vy, p, ax=axes[1, ind], show=False)
        vlt.plot(vz, p, ax=axes[2, ind], show=False)
        vlt.plot(mag_ne, p, ax=axes[3, ind], show=False)

    plt.suptitle(title)
    vlt.auto_adjust_subplots(subplot_params=dict(top=0.9, right=0.9))
    plt.gcf().set_size_inches(6, 7)

    plt.savefig(next_plot_fname(__file__))
    if show:
        vlt.mplshow()

    def test_func(self):

        @cuda.jit(argtypes=[float32[:, ::1], float32[:, ::1], float32[:, ::1]])
        def cu_square_matrix_mul(A, B, C):
            sA = cuda.shared.array(shape=SM_SIZE, dtype=float32)
            sB = cuda.shared.array(shape=(tpb, tpb), dtype=float32)

            tx = cuda.threadIdx.x
            ty = cuda.threadIdx.y
            bx = cuda.blockIdx.x
            by = cuda.blockIdx.y
            bw = cuda.blockDim.x
            bh = cuda.blockDim.y

            x = tx + bx * bw
            y = ty + by * bh

            acc = float32(0)  # forces all the math to be f32
            for i in range(bpg):
                if x < n and y < n:
                    sA[ty, tx] = A[y, tx + i * tpb]
                    sB[ty, tx] = B[ty + i * tpb, x]

                cuda.syncthreads()

                if x < n and y < n:
                    for j in range(tpb):
                        acc += sA[ty, j] * sB[j, tx]

                cuda.syncthreads()

            if x < n and y < n:
                C[y, x] = acc

        np.random.seed(42)
        A = np.array(np.random.random((n, n)), dtype=np.float32)
        B = np.array(np.random.random((n, n)), dtype=np.float32)
        C = np.empty_like(A)

        s = time()
        stream = cuda.stream()
        with stream.auto_synchronize():
            dA = cuda.to_device(A, stream)
            dB = cuda.to_device(B, stream)
            dC = cuda.to_device(C, stream)
            cu_square_matrix_mul[(bpg, bpg), (tpb, tpb), stream](dA, dB, dC)
            dC.copy_to_host(C, stream)

        e = time()
        tcuda = e - s

        # Host compute
        s = time()
        Cans = np.dot(A, B)
        e = time()
        tcpu = e - s

        # Check result
        np.testing.assert_allclose(C, Cans, rtol=1e-5)

def gaussian_original(params, n_rep=DEFAULT_TRIALS, rec_time=False):  
	'''
	Generates n_rep number of Guassian facilitation curves for Go response for all 
	simulated trials required
	
	Parameters
	-------------
	params : sequence (4,) of float
		a_facGo - amplitude of gaussian curve
		b_facGo - time to peak of gaussian curve
		c_facGo - curvature of gaussian curve
		
		Returns
		--------
		fac_i : array
			facilitation curves for all simulated trials
		t : array
			sequence of time index
	'''

	# expand params
	a_facGo_mean, a_facGo_sd, b_facGo_mean, b_facGo_sd, c_facGo_mean, c_facGo_sd, \
	inhib_mean, inhib_sd = params

	t = np.linspace(-.4, .2, 600, endpoint=False)
	#tau_facGo = 2  # Currently set, but will need to optomize

	# generates n_rep random numbers from a normal distribution of mean, sd that given into function
	a_facGo = np.random.normal(a_facGo_mean, a_facGo_sd, size=n_rep)
	b_facGo = np.random.normal(b_facGo_mean, b_facGo_sd, size=n_rep)
	c_facGo = np.random.normal(c_facGo_mean, c_facGo_sd, size=n_rep)


	# had to change from fac_i, t - why does this cause error now?!?! 
	# sets up empty array of zeros for all simulated trials
	fac_i = np.zeros((n_rep, t.size))
 

	inhib_tonic = np.zeros((n_rep, t.size))	
	inhib = np.random.normal(inhib_mean, inhib_sd, size=n_rep)
	inhib_tonic += inhib[:,np.newaxis]
		
	# time performance if required
	if (rec_time):
		t_start = time()

	for i in range(n_rep):  # for each simulated trial
		# takes parameters passed into model plus pre_t number randomly generated
		# for that simulated trial
		myparams_fac = a_facGo[i], b_facGo[i], c_facGo[i]  
		fac_i[i] = _gaussian_original(t, myparams_fac)  # generates curve for that simulated trial
		#inhib_tonic[i] = get_inhib_tonic(t, inhib[i])
	
	if (rec_time):
		t_diff = time() - t_start
		tps = n_rep / t_diff
		print "Original trials per second: %.0f" % (tps)

	return fac_i, inhib_tonic, t

def timeit(f, *args, **kwargs):
    from timeit import default_timer as time
    t0 = time()
    ret = f(*args, **kwargs)
    t1 = time()

    print("Took {0:.03g} secs.".format(t1 - t0))
    return ret

def timereps(reps, func, *args, **kwargs):
    arr = [None] * reps
    for i in range(reps):
        start = time()
        func(*args, **kwargs)
        end = time()
        arr[i] = end - start
    return min(arr), max(arr), sum(arr) / reps