 	double gamma, sigma;
 	double omega;
 	double a, b, c1, c0;
+
+	cudaSetDevice(0);
 
 	MatrixC = m_new(1, 2 * filters);
 	MatrixA = m_new(2 * filters, 2 * filters); //BLOCKMATRIX

	m_print(state);
#endif

	for(i = 0; i < samples;) {
		block_samples = m_multiply(block_CA,state);


	Matrix MatrixCA;
	Matrix statetmp;

	Matrix dMatrixA_pow, dblock_CA, dstate, dblock_samples;

	MatrixCA = m_multiply(MatrixC, MatrixA);

#if DEBUG == 2
