Skip to content

P
pytensor
提交 2413550e authored 作者: James Bergstra's avatar James Bergstra
浏览文件
操作

added longer lenet_32 test

上级 44903030 132c9d49
隐藏空白字符变更
内嵌 并排
正在显示 包含 654 行增加 和 79 行删除
basic_ops.py
浏览文件 @ 2413550e
......@@ -283,13 +283,15 @@ class GpuDimShuffle(Op):
#reassign the dimension and strides in the host pointers
for i, o in enumerate(self.new_order):
if o == 'x':
#TODO: remove this assertion
# the correct thing to do is to insert a run-time check
# that the size in this dimension is 1
assert node.outputs[0].type.broadcastable[i]
print >> sio, """
CudaNdarray_set_dim(cnda_%(res)s, %(i)s, 1);
CudaNdarray_set_stride(cnda_%(res)s, %(i)s, 0);
""" %locals()
else:
assert not node.outputs[0].type.broadcastable[i]
print >> sio, """
CudaNdarray_set_dim(cnda_%(res)s, %(i)s, CudaNdarray_HOST_DIMS(cnda_%(input)s)[%(o)s]);
CudaNdarray_set_stride(cnda_%(res)s, %(i)s, CudaNdarray_HOST_STRIDES(cnda_%(input)s)[%(o)s]);
......@@ -335,6 +337,26 @@ class GpuDimShuffle(Op):
return (1,0)
class GpuSum(Op):
"""GpuSum is a Reduction along some dimensions by summation.
The dimensions along which to sum is specified by the `reduce_mask` that you pass to the
constructor. The `reduce_mask` is a tuple of booleans (actually integers 0 or 1) that
specify for each input dimension, whether to reduce it (1) or not (0).
For example:
- reduce_mask == (1,) sums a vector to a scalar
- reduce_mask == (1,0) computes the sum of each column in a matrix
- reduce_mask == (0,1) computes the sum of each row in a matrix
- reduce_mask == (1,1,1) computes the sum of all elements in a 3-tensor.
:note: any reduce_mask of all zeros is a sort of 'copy', and may be removed during graph
optimization
"""
def __init__(self, reduce_mask):
self.reduce_mask = tuple(reduce_mask)
......@@ -356,6 +378,435 @@ class GpuSum(Op):
def perform(self, node, (x,), (z,)):
z[0] = x.reduce_sum(self.reduce_mask)
def c_code(self, node, name, (x,), (z,), sub):
nd_in = node.inputs[0].type.ndim
nd_out = node.outputs[0].type.ndim
assert nd_in - nd_out == sum(self.reduce_mask)
sio = StringIO.StringIO()
fail = sub['fail']
#check input
print >> sio, """
if (cnda_%(x)s->nd != %(nd_in)s)
{
PyErr_Format(PyExc_TypeError, "required nd=%(nd_in)s, got nd=%%i", cnda_%(x)s->nd);
%(fail)s;
}
""" %locals()
#
# alloc an output if we need one
#
# check the basics of out output
print >> sio, """
if ( !cnda_%(z)s
|| (cnda_%(z)s->nd != %(nd_out)s)
""" % locals()
#ensure that the output has the right non-reduced dimensions
j = 0
for i in xrange(nd_in):
if not self.reduce_mask[i]:
print >> sio, " || (CudaNdarray_HOST_DIMS(cnda_%(z)s)[%(j)s] !=CudaNdarray_HOST_DIMS(cnda_%(x)s)[%(i)s]) " % locals()
j += 1
print >> sio, """
)
{
""" %locals()
print >> sio, "int new_dims[%(nd_out)s]; " % locals()
j = 0
for i in xrange(nd_in):
if not self.reduce_mask[i]:
print >> sio, 'new_dims[%(j)s] = CudaNdarray_HOST_DIMS(cnda_%(x)s)[%(i)s];' % locals()
j += 1
print >> sio, """
Py_XDECREF(cnda_%(z)s);
cnda_%(z)s = (CudaNdarray*) CudaNdarray_NewDims(%(nd_out)s, new_dims);
if (NULL == cnda_%(z)s)
{
%(fail)s;
}
}
""" %locals()
#
# Now perform the reduction
#
if self.reduce_mask == (1,):
self.c_code_reduce_1(sio, node, name, x, z, fail)
elif self.reduce_mask == (1,1):
self.c_code_reduce_11(sio, node, name, x, z, fail)
elif self.reduce_mask == (1,0):
self.c_code_reduce_10(sio, node, name, x, z, fail)
elif self.reduce_mask == (1,0,1,1):
self.c_code_reduce_1011(sio, node, name, x, z, fail)
else:
print 'UNWRITTEN REDUCE MASK', self.reduce_mask
assert 0
return sio.getvalue()
def c_code_reduce_1(self, sio, node, name, x, z, fail):
print >> sio, """
{
int verbose = 0;
dim3 n_threads(
std::min(CudaNdarray_HOST_DIMS(cnda_%(x)s)[0],
NUM_VECTOR_OP_THREADS_PER_BLOCK));
dim3 n_blocks(1);
if (verbose) printf("running kernel_reduce_sum_1_%(name)s\\n");
int n_shared = sizeof(float) * n_threads.x * n_threads.y * n_threads.z;
kernel_reduce_sum_1_%(name)s<<<n_blocks, n_threads, n_shared>>>(
CudaNdarray_HOST_DIMS(cnda_%(x)s)[0],
CudaNdarray_DEV_DATA(cnda_%(x)s),
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[0],
CudaNdarray_DEV_DATA(cnda_%(z)s));
CNDA_THREAD_SYNC;
if (cudaSuccess != cudaGetLastError())
{
%(fail)s;
}
}
""" %locals()
def c_code_reduce_11(self, sio, node, name, x, z, fail):
print >> sio, """
{
int verbose = 0;
dim3 n_threads(
std::min(CudaNdarray_HOST_DIMS(cnda_%(x)s)[1],
NUM_VECTOR_OP_THREADS_PER_BLOCK));
while (n_threads.y * n_threads.x < NUM_VECTOR_OP_THREADS_PER_BLOCK) ++n_threads.y;
n_threads.y -= 1;
if (n_threads.y > CudaNdarray_HOST_DIMS(cnda_%(x)s)[0])
n_threads.y = CudaNdarray_HOST_DIMS(cnda_%(x)s)[0];
dim3 n_blocks(1);
if (verbose) printf("running kernel_reduce_sum_11_%(name)s\\n");
int n_shared = sizeof(float) * n_threads.x * n_threads.y * n_threads.z;
kernel_reduce_sum_11_%(name)s<<<n_blocks, n_threads, n_shared>>>(
CudaNdarray_HOST_DIMS(cnda_%(x)s)[0],
CudaNdarray_HOST_DIMS(cnda_%(x)s)[1],
CudaNdarray_DEV_DATA(cnda_%(x)s),
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[0],
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[1],
CudaNdarray_DEV_DATA(cnda_%(z)s));
CNDA_THREAD_SYNC;
if (cudaSuccess != cudaGetLastError())
{
%(fail)s;
}
}
""" %locals()
def c_code_reduce_10(self, sio, node, name, x, z, fail):
print >> sio, """
{
int verbose = 0;
dim3 n_threads(
std::min(CudaNdarray_HOST_DIMS(cnda_%(x)s)[0],
NUM_VECTOR_OP_THREADS_PER_BLOCK));
dim3 n_blocks(CudaNdarray_HOST_DIMS(cnda_%(x)s)[1]);
if (verbose) printf("running kernel_reduce_sum_10_%(name)s\\n");
int n_shared = sizeof(float) * n_threads.x;
kernel_reduce_sum_10_%(name)s<<<n_blocks, n_threads, n_shared>>>(
CudaNdarray_HOST_DIMS(cnda_%(x)s)[0],
CudaNdarray_HOST_DIMS(cnda_%(x)s)[1],
CudaNdarray_DEV_DATA(cnda_%(x)s),
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[0],
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[1],
CudaNdarray_DEV_DATA(cnda_%(z)s),
CudaNdarray_HOST_STRIDES(cnda_%(z)s)[0]
);
CNDA_THREAD_SYNC;
if (cudaSuccess != cudaGetLastError())
{
%(fail)s;
}
}
""" %locals()
def c_code_reduce_1011(self, sio, node, name, x, z, fail):
print >> sio, """
{
int verbose = 0;
dim3 n_threads(
std::min(CudaNdarray_HOST_DIMS(cnda_%(x)s)[3],
NUM_VECTOR_OP_THREADS_PER_BLOCK));
while (n_threads.y * n_threads.x < NUM_VECTOR_OP_THREADS_PER_BLOCK) ++n_threads.y;
n_threads.y -= 1;
if (n_threads.y > CudaNdarray_HOST_DIMS(cnda_%(x)s)[2])
n_threads.y = CudaNdarray_HOST_DIMS(cnda_%(x)s)[2];
while (n_threads.x * n_threads.y * n_threads.z < NUM_VECTOR_OP_THREADS_PER_BLOCK) ++n_threads.z;
n_threads.z -= 1;
if (n_threads.z > 64)
n_threads.z = 64;
if (n_threads.z > CudaNdarray_HOST_DIMS(cnda_%(x)s)[0])
n_threads.z = CudaNdarray_HOST_DIMS(cnda_%(x)s)[0];
dim3 n_blocks(CudaNdarray_HOST_DIMS(cnda_%(x)s)[1]);
if (verbose) printf("running kernel_reduce_sum_1011_%(name)s\\n");
if (verbose) fprint_CudaNdarray(stdout, cnda_%(x)s);
if (verbose) fprint_CudaNdarray(stdout, cnda_%(z)s);
int n_shared = sizeof(float) * n_threads.x * n_threads.y * n_threads.z;
kernel_reduce_sum_1011_%(name)s<<<n_blocks, n_threads, n_shared>>>(
CudaNdarray_HOST_DIMS(cnda_%(x)s)[0],
CudaNdarray_HOST_DIMS(cnda_%(x)s)[1],
CudaNdarray_HOST_DIMS(cnda_%(x)s)[2],
CudaNdarray_HOST_DIMS(cnda_%(x)s)[3],
CudaNdarray_DEV_DATA(cnda_%(x)s),
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[0],
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[1],
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[2],
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[3],
CudaNdarray_DEV_DATA(cnda_%(z)s),
CudaNdarray_HOST_STRIDES(cnda_%(z)s)[0]);
CNDA_THREAD_SYNC;
if (cudaSuccess != cudaGetLastError())
{
%(fail)s;
}
}
""" %locals()
def c_code_cache_version(self):
return ()
def c_support_code_apply(self, node, nodename):
sio = StringIO.StringIO()
if self.reduce_mask == (1,):
#this kernel is ok for up to a few thousand elements, but
# it only runs on ONE multiprocessor
print >> sio, """
static __global__ void kernel_reduce_sum_1_%(nodename)s(
const unsigned int d0,
const float *A, const int sA0,
float * Z)
{
const int threadCount = blockDim.x;
const int threadNum = threadIdx.x;
extern __shared__ float buf[];
float mysum = 0.0f;
if (warpSize != 32)
{
return; //TODO: set error code
}
for (int i0 = threadIdx.x; i0 < d0; i0 += blockDim.x)
{
float Ai = A[i0 * sA0];
mysum += Ai;
}
buf[threadNum] = mysum;
__syncthreads();
// rest of function is handled by one warp
if (threadNum < warpSize)
{
for (int i = threadNum + warpSize; i < threadCount; i += warpSize)
{
mysum += buf[i];
}
buf[threadNum] = mysum;
if (threadNum < 16)
{
//reduce so that threadNum 0 has the sum of everything
if(threadNum + 16 < threadCount) buf[threadNum] += buf[threadNum+16];
if(threadNum + 8 < threadCount) buf[threadNum] += buf[threadNum+8];
if(threadNum + 4 < threadCount) buf[threadNum] += buf[threadNum+4];
if(threadNum + 2 < threadCount) buf[threadNum] += buf[threadNum+2];
if(threadNum + 1 < threadCount) buf[threadNum] += buf[threadNum+1];
if (threadNum == 0)
{
Z[0] = buf[0];
}
}
}
}
""" %locals()
if self.reduce_mask == (1,1):
#this kernel is ok for up to a few thousand elements, but
# it only runs on ONE multiprocessor
print >> sio, """
static __global__ void kernel_reduce_sum_11_%(nodename)s(
const int d0,
const int d1,
const float *A, const int sA0, const int sA1,
float * Z)
{
const int threadCount = blockDim.x * blockDim.y;
const int threadNum = threadIdx.y*blockDim.x + threadIdx.x;
extern __shared__ float buf[];
float mysum = 0.0f;
if (warpSize != 32)
{
return; //TODO: set error code
}
for (int i0 = threadIdx.y; i0 < d0; i0 += blockDim.y)
{
for (int i1 = threadIdx.x; i1 < d1; i1 += blockDim.x)
{
float Ai = A[i0 * sA0 + i1 * sA1];
mysum += Ai;
}
}
buf[threadNum] = mysum;
__syncthreads();
// rest of function is handled by one warp
if (threadNum < warpSize)
{
for (int i = threadNum + warpSize; i < threadCount; i += warpSize)
{
mysum += buf[i];
}
buf[threadNum] = mysum;
if (threadNum < 16)
{
//reduce so that threadNum 0 has the sum of everything
if(threadNum + 16 < threadCount) buf[threadNum] += buf[threadNum+16];
if(threadNum + 8 < threadCount) buf[threadNum] += buf[threadNum+8];
if(threadNum + 4 < threadCount) buf[threadNum] += buf[threadNum+4];
if(threadNum + 2 < threadCount) buf[threadNum] += buf[threadNum+2];
if(threadNum + 1 < threadCount) buf[threadNum] += buf[threadNum+1];
if (threadNum == 0)
{
Z[0] = buf[0];
}
}
}
}
""" %locals()
if self.reduce_mask == (1,0):
# this kernel uses one block for each column,
# threads per block for each element per column.
#TODO: This kernel is pretty inefficient in terms of reading, because if A is
# c_contiguous (typical case) then each warp is accessing non-contigous
# memory (a segment of a column).
print >> sio, """
static __global__ void kernel_reduce_sum_10_%(nodename)s(
const int d0,
const int d1,
const float *A, const int sA0, const int sA1,
float * Z, const int sZ0)
{
const int threadCount = blockDim.x;
const int threadNum = threadIdx.x;
extern __shared__ float buf[];
float mysum = 0.0f;
if (warpSize != 32)
{
return; //TODO: set error code
}
for (int i0 = threadIdx.x; i0 < d0; i0 += blockDim.x)
{
float Ai = A[i0 * sA0 + blockIdx.x * sA1];
mysum += Ai;
}
buf[threadNum] = mysum;
__syncthreads();
// rest of function is handled by one warp
if (threadNum < warpSize)
{
for (int i = threadNum + warpSize; i < threadCount; i += warpSize)
{
mysum += buf[i];
}
buf[threadNum] = mysum;
if (threadNum < 16)
{
//reduce so that threadNum 0 has the sum of everything
if(threadNum + 16 < threadCount) buf[threadNum] += buf[threadNum+16];
if(threadNum + 8 < threadCount) buf[threadNum] += buf[threadNum+8];
if(threadNum + 4 < threadCount) buf[threadNum] += buf[threadNum+4];
if(threadNum + 2 < threadCount) buf[threadNum] += buf[threadNum+2];
if(threadNum + 1 < threadCount) buf[threadNum] += buf[threadNum+1];
if (threadNum == 0)
{
Z[blockIdx.x * sZ0] = buf[0];
}
}
}
}
""" %locals()
if self.reduce_mask == (1,0,1,1):
print >> sio, """
static __global__ void kernel_reduce_sum_1011_%(nodename)s(
const unsigned int d0,
const unsigned int d1,
const unsigned int d2,
const unsigned int d3,
const float *A, const int sA0, const int sA1, const int sA2, const int sA3,
float * Z, const int sZ0)
{
const int threadCount = blockDim.x * blockDim.y * blockDim.z;
const int threadNum = threadIdx.z * blockDim.x * blockDim.y + threadIdx.y * blockDim.x + threadIdx.x;
extern __shared__ float buf[];
float mysum = 0.0f;
if (warpSize != 32)
{
return; //TODO: set error code
}
for (int i0 = threadIdx.z; i0 < d0; i0 += blockDim.z)
{
for (int i2 = threadIdx.y; i2 < d2; i2 += blockDim.y)
{
for (int i3 = threadIdx.x; i3 < d3; i3 += blockDim.x)
{
float Ai = A[i0 * sA0 + blockIdx.x * sA1 + i2 * sA2 + i3 * sA3];
mysum += Ai;
}
}
}
buf[threadNum] = mysum;
__syncthreads();
// rest of function is handled by one warp
if (threadNum < warpSize)
{
for (int i = threadNum + warpSize; i < threadCount; i += warpSize)
{
mysum += buf[i];
}
buf[threadNum] = mysum;
if (threadNum < 16)
{
//reduce so that threadNum 0 has the sum of everything
if(threadNum + 16 < threadCount) buf[threadNum] += buf[threadNum+16];
if(threadNum + 8 < threadCount) buf[threadNum] += buf[threadNum+8];
if(threadNum + 4 < threadCount) buf[threadNum] += buf[threadNum+4];
if(threadNum + 2 < threadCount) buf[threadNum] += buf[threadNum+2];
if(threadNum + 1 < threadCount) buf[threadNum] += buf[threadNum+1];
if (threadNum == 0)
{
Z[blockIdx.x*sZ0] = buf[0];
}
}
}
}
""" %locals()
return sio.getvalue()
class GpuReshape(tensor.Reshape):
# __hash__, __eq__, __str__ come from tensor.Subtensor
def make_node(self, x, shp):
......@@ -375,6 +826,25 @@ class GpuSubtensor(tensor.Subtensor):
return rval
def perform(self, node, inputs, (out, )):
x = inputs[0]
indices = list(reversed(inputs[1:]))
def convert(entry):
if isinstance(entry, Type):
return indices.pop()
elif isinstance(entry, slice):
return slice(convert(entry.start),
convert(entry.stop),
convert(entry.step))
else:
return entry
cdata = tuple(map(convert, self.idx_list))
if len(cdata) == 1:
cdata = cdata[0]
out[0] = x.__getitem__(cdata)
def old_perform(self, node, inputs, (out, )):
indices = list(reversed(inputs[1:]))
def convert(entry):
......
blas.py
浏览文件 @ 2413550e
......@@ -174,9 +174,11 @@ class GpuConv(Op):
def make_node(self, img, kern):
if img.type.ndim != 4:
raise TypeError('img must be 4D tensor')
if img.type != kern.type:
raise TypeError('img and kern must have same type')
return Apply(self, [img, kern], [img.type()])
if kern.type.ndim != 4:
raise TypeError('kern must be 4D tensor')
broadcastable = [img.type.broadcastable[0], kern.type.broadcastable[0], False, False]
return Apply(self, [img, kern], [CudaNdarrayType(broadcastable)()])
def perform(self, node, (img, kern), (out,)):
out[0] = cuda_ndarray.conv(img, kern,
......@@ -187,13 +189,28 @@ class GpuConv(Op):
kern_align=self.logical_kern_align_top,
verbose=0)
from theano.sandbox.downsample import DownsampleFactorMax
class GpuDownsampleFactorMax(DownsampleFactorMax):
# inherit __eq__, __hash__, __str__
class GpuDownsampleFactorMax(Op):
def __init__(self, ds, ignore_border=False):
self.ds = tuple(ds)
self.ignore_border = ignore_border
def __eq__(self, other):
return type(self) == type(other) and self.ds == other.ds and self.ignore_border == other.ignore_border
def __hash__(self):
return hash(type(self)) ^ hash(self.ds) ^ hash(self.ignore_border)
def __str__(self):
return '%s{%s,%s}' % (self.__class__.__name__, self.ds, self.ignore_border)
def make_node(self, x):
if not isinstance(x.type, CudaNdarrayType):
raise TypeError()
if not x.type.ndim == 4:
raise TypeError()
return Apply(self, [x], [x.type()])
def perform(self, node, input_storage, output_storage):
raise NotImplementedError('only C is implemented')
#def perform(self, node, input_storage, output_storage):
#raise NotImplementedError('only C is implemented')
def c_code_cache_version(self):
return ()
def c_code(self, node, nodename, (x,), (z,), sub):
......@@ -240,8 +257,8 @@ class GpuDownsampleFactorMax(DownsampleFactorMax):
//dim3 block(std::min(dims[3], 512)); //TODO: implement this by supporting more
//outputs than threads
dim3 block(dims[3]);
int shared= xdim3*sizeof(float);
kMaxPool_%(nodename)s<%(ds0)s, %(ds1)s> <<<grid, block, shared>>>(
if ((grid.x*grid.y) && dims[3])
kMaxPool_%(nodename)s<%(ds0)s, %(ds1)s> <<<grid, block, xdim3*sizeof(float)>>>(
dims[0], dims[1], dims[2], dims[3], xdim2, xdim3,
CudaNdarray_DEV_DATA(cnda_%(x)s),
CudaNdarray_HOST_STRIDES(cnda_%(x)s)[0],
......@@ -253,8 +270,14 @@ class GpuDownsampleFactorMax(DownsampleFactorMax):
cudaError_t err = cudaGetLastError();
if( cudaSuccess != err)
{
PyErr_Format(PyExc_RuntimeError, "Cuda error: %%s: %%s.threads.x=%%d threads.y=%%d threads.z=%%d grid.x=%%d grid.y=%%d shared=%%d\\n", "kMaxPool_%(nodename)s",
cudaGetErrorString(err), block.x, block.y, block.z, grid.x, grid.y, shared);
PyErr_Format(PyExc_RuntimeError, "Cuda error: %%s: %%s. (grid: %%i x %%i; block: %%i x %%i x %%i)\\n",
"kMaxPool_%(nodename)s",
cudaGetErrorString(err),
grid.x,
grid.y,
block.x,
block.y,
block.z);
%(fail)s;
}
}
......@@ -270,8 +293,8 @@ class GpuDownsampleFactorMax(DownsampleFactorMax):
float *z)
{
float cur_max, cur_x;
int i0 = blockIdx.x / D0;
int i1 = blockIdx.x %% D0;
int i0 = blockIdx.x %% D0;
int i1 = blockIdx.x / D0;
int i2 = blockIdx.y;
extern __shared__ float xbuf[]; //size [xD3]
......@@ -280,9 +303,9 @@ class GpuDownsampleFactorMax(DownsampleFactorMax):
{
__syncthreads();
// load the current row of the image into shared memory
for (int i3 = threadIdx.x; i3 < xD3; i3 += blockDim.x)
for (int j = threadIdx.x; j < xD3; j += blockDim.x)
{
xbuf[i3] = x[i0*xS0 + i1*xS1 + (i2*pf2+r2)*xS2 + i3*xS3];
xbuf[j] = x[i0*xS0 + i1*xS1 + (i2*pf2+r2)*xS2 + j*xS3];
}
__syncthreads();
......@@ -290,10 +313,24 @@ class GpuDownsampleFactorMax(DownsampleFactorMax):
cur_max = (r2 == 0) ? xbuf[threadIdx.x*pf3] : cur_max;
// do a mini-reduction over the pf3 relevant elements in the current row
for (int k = 0; k < pf3; ++k)
if (%(ignore_border)s)
{
for (int k = 0; k < pf3; ++k)
{
cur_x = xbuf[threadIdx.x*pf3+k];
cur_max = (cur_x > cur_max) ? cur_x : cur_max;
}
}
else
{
cur_x = xbuf[threadIdx.x*pf3+k];
cur_max = (cur_x < cur_max) ? cur_x : cur_max;
for (int k = 0; k < pf3; ++k)
{
if (threadIdx.x*pf3 + k < xD3)
{
cur_x = xbuf[threadIdx.x*pf3+k];
cur_max = (cur_x > cur_max) ? cur_x : cur_max;
}
}
}
}
......@@ -302,13 +339,24 @@ class GpuDownsampleFactorMax(DownsampleFactorMax):
}
""" % locals()
from theano.sandbox.downsample import DownsampleFactorMaxGrad
class GpuDownsampleFactorMaxGrad(DownsampleFactorMaxGrad):
# inherit __eq__, __hash__, __str__
class GpuDownsampleFactorMaxGrad(Op):
def __init__(self, ds, ignore_border):
self.ds = tuple(ds)
self.ignore_border = ignore_border
def __eq__(self, other):
return type(self) == type(other) and self.ds == other.ds and self.ignore_border == other.ignore_border
def __hash__(self):
return hash(type(self)) ^ hash(self.ds) ^ hash(self.ignore_border)
def __str__(self):
return '%s{%s,%s}' % (self.__class__.__name__, self.ds, self.ignore_border)
def make_node(self, x, z, gz):
return Apply(self, [x, z, gz], [x.type()])
def perform(self, node, input_storage, output_storage):
raise NotImplementedError('only C is implemented')
#def perform(self, node, input_storage, output_storage):
#raise NotImplementedError('only C is implemented')
def c_code_cache_version(self):
return ()
def c_code(self, node, nodename, (x, z, gz), (gx,), sub):
......@@ -340,9 +388,9 @@ class GpuDownsampleFactorMaxGrad(DownsampleFactorMaxGrad):
}
}
{
dim3 grid(CudaNdarray_HOST_DIMS(cnda_%(x)s)[0], CudaNdarray_HOST_DIMS(cnda_%(x)s)[2]);
//TODO: implement this by supporting more
//outputs than threads
dim3 grid(CudaNdarray_HOST_DIMS(cnda_%(x)s)[0], CudaNdarray_HOST_DIMS(cnda_%(x)s)[2]);
dim3 block(CudaNdarray_HOST_DIMS(cnda_%(x)s)[3]);
kDownsampleMaxGrad_%(nodename)s<%(ds0)s, %(ds1)s> <<<grid, block>>>(
CudaNdarray_HOST_DIMS(cnda_%(z)s)[0],
......@@ -401,9 +449,11 @@ class GpuDownsampleFactorMaxGrad(DownsampleFactorMaxGrad):
int i2 = blockIdx.y; // row wrt z and/or gz
int x_col = threadIdx.x;
// The algorithm here is that every thread writes one output pixel per line
//TODO: raise occupancy. Use threadIdx.y to run several iterations of this i1 loop
//in parallel
for (i1 = 0; i1 < D1; ++i1)
{
// The algorithm here is that every thread writes one output pixel per line
if (%(ignore_border)s && (x_col >= ds1 * D3))
{
my_gz = 0;
......@@ -417,7 +467,7 @@ class GpuDownsampleFactorMaxGrad(DownsampleFactorMaxGrad):
for (int x_row = i2*ds0; (x_row < i2*ds0+ds0) && (%(ignore_border)s || (x_row < xD2)); ++x_row)
{
gx[i0 * D1*xD2*xD3 + i1*xD2*xD3 + x_row*xD3 + x_col]
= (my_z == x[i0*xS0 + i1*xS1 + x_row*xS2 + x_col]) ? my_gz : 0;
= (my_z == x[i0*xS0 + i1*xS1 + x_row*xS2 + x_col*xS3]) ? my_gz : 0;
}
}
}
......
nnet.py
浏览文件 @ 2413550e
......@@ -186,9 +186,6 @@ class GpuCrossentropySoftmax1HotWithBiasDx (Op):
return self.__class__.__name__
def make_node(self, dy, sm, y_idx):
return Apply(self, [dy, sm, y_idx],[sm.type()])
def perform(self, node, input_storage, output_storage):
raise NotImplementedError('only C is implemented')
def c_code_cache_version(self):
return ()
def c_code(self, node, nodename, (dnll, sm, y_idx), (dx,), sub):
......
tests/test_blas.py
浏览文件 @ 2413550e
......@@ -7,6 +7,7 @@ import numpy
import theano_cuda_ndarray as tcn
from theano.sandbox.downsample import DownsampleFactorMax
def test_dot():
......@@ -46,20 +47,66 @@ def test_gemm():
assert numpy.allclose(numpy.dot(a0, bval)+numpy.exp(cval), a.value)
def test_maxpool():
"""TODO: test the gpu version!!! """
for d0, d1, r_true, r_false in [(4,4,[[[[5,7],[13,15]]]],[[[[5,7],[13,15]]]]),
(5,5,[[[[6, 8],[ 16, 18], [ 21, 23]]]],
[[[[6, 8, 9],[ 16, 18, 19], [ 21, 23, 24]]]])]:
for border,ret in [(True,r_true),(False, r_false)]:
ret=numpy.array(ret)
a = tcn.blas.DownsampleFactorMax((2,2),border)
dmatrix4 = tensor.TensorType("float32", (False, False, False, False))
b = dmatrix4()
f = pfunc([b], [a(b)])
bval = numpy.arange(0,d0*d1).reshape(1,1,d0,d1)
r = f(bval)[0]
# print bval, bval.shape, border
print r, r.shape
assert (ret==r).all()
if 0:
# This is commented out because it doesn't make sense...
# tcn.blas has no op called DownsampleFactorMax
# tcn.blas has an op called GpuDownsampleFactorMax, but that op requires arguments that are
# CudaNdarrayType variables... so rethink this test?
def test_maxpool():
"""TODO: test the gpu version!!! """
for d0, d1, r_true, r_false in [(4,4,[[[[5,7],[13,15]]]],[[[[5,7],[13,15]]]]),
(5,5,[[[[6, 8],[ 16, 18], [ 21, 23]]]],
[[[[6, 8, 9],[ 16, 18, 19], [ 21, 23, 24]]]])]:
for border,ret in [(True,r_true),(False, r_false)]:
ret=numpy.array(ret)
a = tcn.blas.DownsampleFactorMax((2,2),border)
dmatrix4 = tensor.TensorType("float32", (False, False, False, False))
b = dmatrix4()
f = pfunc([b], [a(b)])
bval = numpy.arange(0,d0*d1).reshape(1,1,d0,d1)
r = f(bval)[0]
# print bval, bval.shape, border
print r, r.shape
assert (ret==r).all()
def test_downsample():
for shp in [
(1, 1, 1, 12),
(1, 1, 2, 2),
#(1, 1, 1, 1), #### Commented out because it makes FP-exception that I don't understand
(1,1,4,4),
(1, 1, 10, 11),
(1, 2, 2, 2),
(3,5,4,4),
(1, 1, 12, 12),
(1, 1, 2, 14),
(1, 1, 12, 14),
(1, 1, 14, 14),
(1, 1, 16, 16),
(1, 1, 18, 18),
(1, 1, 24, 24),
(1, 6, 24, 24),
(10, 1, 24, 24),
(10, 6, 24, 24),
(30, 6, 12, 12),
(30, 2, 24, 24),
(30, 6, 24, 24),
(10, 10, 10, 11)]:
for ds in (1,1), (2, 2):
if ds[0] > shp[2]: continue
if ds[1] > shp[3]: continue
for ignore_border in (True, False):
print 'test_downsample', shp, ds, ignore_border
ds_op = DownsampleFactorMax(ds, ignore_border=ignore_border)
a = tcn.shared_constructor(numpy.random.rand(*shp), 'a')
f = pfunc([], ds_op(tensor.as_tensor_variable(a)))
worked = False
for i, node in enumerate(f.maker.env.toposort()):
print i, node
if isinstance(node.op, tcn.blas.GpuDownsampleFactorMax):
f() # let debugmode do the testing
worked = True
assert worked
tests/test_nnet.py
浏览文件 @ 2413550e
......@@ -14,7 +14,7 @@ import numpy
import theano_cuda_ndarray as tcn
import logging
logging.getLogger('theano.gradient').setLevel(logging.INFO)
logging.getLogger('test_cuda_ndarray.tests.test_nnet').setLevel(logging.INFO)
def get_mode():
......@@ -97,18 +97,18 @@ def run_conv_nnet1(shared_fn):
n_out = 10
w = shared_fn(numpy.asarray(0.01*(numpy.random.rand(*shape_kern)-0.5), dtype='float32'), 'w')
b = shared_fn(numpy.asarray(numpy.zeros((n_kern,1,1)), dtype='float32'), 'b')
b = shared_fn(numpy.asarray(numpy.zeros((n_kern,)), dtype='float32'), 'b')
v = shared_fn(numpy.asarray(numpy.zeros((n_hid, n_out)), dtype='float32'), 'c')
c = shared_fn(numpy.asarray(numpy.zeros(n_out), dtype='float32'), 'c')
x = tensor.Tensor(dtype='float32', broadcastable=(0,0,0,0))('x')
x = tensor.Tensor(dtype='float32', broadcastable=(0,1,0,0))('x')
y = tensor.fmatrix('y')
lr = tensor.fscalar('lr')
conv_op = theano.sandbox.conv.ConvOp(shape_img[2:], shape_kern[2:], n_kern, n_batch, 1, 1)
conv_op.set_flops()
hid = tensor.tanh(conv_op(x, w)+b)
hid = tensor.tanh(conv_op(x, w)+b.dimshuffle((0,'x','x')))
hid_flat = hid.reshape((n_batch, n_hid))
out = tensor.tanh(tensor.dot(hid_flat, v)+c)
loss = tensor.sum(0.5 * (out-y)**2 * lr)
......@@ -174,13 +174,13 @@ def run_conv_nnet2(shared_fn): # pretend we are training LeNet for MNIST
n_out = 10
w0 = shared_fn(numpy.asarray(0.01*(numpy.random.rand(*shape_kern)-0.5), dtype='float32'), 'w0')
b0 = shared_fn(numpy.asarray(numpy.zeros((n_kern,1,1)), dtype='float32'), 'b0')
b0 = shared_fn(numpy.asarray(numpy.zeros((n_kern,)), dtype='float32'), 'b0')
w1 = shared_fn(numpy.asarray(0.01*(numpy.random.rand(*shape_kern1)-0.5), dtype='float32'), 'w1')
b1 = shared_fn(numpy.asarray(numpy.zeros((n_kern1,1,1)), dtype='float32'), 'b1')
b1 = shared_fn(numpy.asarray(numpy.zeros((n_kern1,)), dtype='float32'), 'b1')
v = shared_fn(numpy.asarray(numpy.zeros((n_hid, n_out)), dtype='float32'), 'c')
c = shared_fn(numpy.asarray(numpy.zeros(n_out), dtype='float32'), 'c')
x = tensor.Tensor(dtype='float32', broadcastable=(0,0,0,0))('x')
x = tensor.Tensor(dtype='float32', broadcastable=(0,1,0,0))('x')
y = tensor.fmatrix('y')
lr = tensor.fscalar('lr')
......@@ -188,10 +188,9 @@ def run_conv_nnet2(shared_fn): # pretend we are training LeNet for MNIST
conv_op1 = theano.sandbox.conv.ConvOp((n_kern,logical_hid_shape[0]/2, logical_hid_shape[1]/2), shape_kern1[2:], n_kern1, n_batch, 1, 1)
conv_op.set_flops()
conv_op1.set_flops()
hid = tensor.tanh(conv_op(x, w0)+b0)
hid1 = tensor.tanh(conv_op1(hid[:,:,::2,::2], w1) + b1)
hid = tensor.tanh(conv_op(x, w0)+b0.dimshuffle((0,'x','x')))
hid1 = tensor.tanh(conv_op1(hid[:,:,::2,::2], w1) + b1.dimshuffle((0,'x','x')))
hid_flat = hid1.reshape((n_batch, n_hid))
out = tensor.tanh(tensor.dot(hid_flat, v)+c)
loss = tensor.sum(0.5 * (out-y)**2 * lr)
......@@ -226,7 +225,7 @@ def test_conv_nnet2():
print rval_cpu[0], rval_gpu[0],rval_cpu[0]-rval_gpu[0]
assert numpy.allclose(rval_cpu, rval_gpu,rtol=1e-4,atol=1e-4)
def run_conv_nnet2_classif(shared_fn, isize, ksize, n_batch=60, n_iter=25):
def run_conv_nnet2_classif(shared_fn, isize, ksize, n_batch, n_iter):
shape_img = (n_batch, 1, isize, isize)
......@@ -243,13 +242,13 @@ def run_conv_nnet2_classif(shared_fn, isize, ksize, n_batch=60, n_iter=25):
n_out = 10
w0 = shared_fn(numpy.asarray(0.01*(numpy.random.rand(*shape_kern)-0.5), dtype='float32'), 'w0')
b0 = shared_fn(numpy.asarray(numpy.zeros((n_kern,1,1)), dtype='float32'), 'b0')
b0 = shared_fn(numpy.asarray(numpy.zeros((n_kern,)), dtype='float32'), 'b0')
w1 = shared_fn(numpy.asarray(0.01*(numpy.random.rand(*shape_kern1)-0.5), dtype='float32'), 'w1')
b1 = shared_fn(numpy.asarray(numpy.zeros((n_kern1,1,1)), dtype='float32'), 'b1')
b1 = shared_fn(numpy.asarray(numpy.zeros((n_kern1,)), dtype='float32'), 'b1')
v = shared_fn(numpy.asarray(0.01*numpy.random.randn(n_hid, n_out), dtype='float32'), 'c')
c = shared_fn(numpy.asarray(numpy.zeros(n_out), dtype='float32'), 'c')
x = tensor.Tensor(dtype='float32', broadcastable=(0,0,0,0))('x')
x = tensor.Tensor(dtype='float32', broadcastable=(0,1,0,0))('x')
y = tensor.fmatrix('y')
lr = tensor.fscalar('lr')
......@@ -260,15 +259,15 @@ def run_conv_nnet2_classif(shared_fn, isize, ksize, n_batch=60, n_iter=25):
ds_op = theano.sandbox.downsample.DownsampleFactorMax((2,2), ignore_border=False)
hid = tensor.tanh(ds_op(conv_op(x, w0)+b0))
hid1 = tensor.tanh(conv_op1(hid, w1) + b1)
hid = tensor.tanh(ds_op(conv_op(x, w0)+b0.dimshuffle((0,'x','x'))))
hid1 = tensor.tanh(conv_op1(hid, w1) + b1.dimshuffle((0,'x','x')))
hid_flat = hid1.reshape((n_batch, n_hid))
out = tensor.nnet.softmax(tensor.dot(hid_flat, v)+c)
loss = tensor.sum(tensor.nnet.crossentropy_categorical_1hot(out, tensor.argmax(y, axis=1)) * lr)
print 'loss type', loss.type
params = [w0, b0, w1, b1, v, c]
gparams = tensor.grad(loss, params)
gparams = tensor.grad(loss, params, warn_type=True)
mode = get_mode()
......@@ -291,16 +290,19 @@ def run_conv_nnet2_classif(shared_fn, isize, ksize, n_batch=60, n_iter=25):
print_mode(mode)
return rvals, t1-t0
def run_test_conv_nnet2_classif(seed, isize, ksize, bsize, ignore_error=False, gpu_only=False):
def cmp_run_conv_nnet2_classif(seed, isize, ksize, bsize,
ignore_error=False,
n_iter=10,
gpu_only=False):
if gpu_only:
numpy.random.seed(seed)
rval_gpu, t = run_conv_nnet2_classif(tcn.shared_constructor, isize, ksize, bsize)
return
numpy.random.seed(seed)
rval_cpu, tc = run_conv_nnet2_classif(shared, isize, ksize, bsize)
rval_gpu, tg = run_conv_nnet2_classif(tcn.shared_constructor, isize, ksize, bsize, n_iter)
numpy.random.seed(seed)
rval_gpu, tg = run_conv_nnet2_classif(tcn.shared_constructor, isize, ksize, bsize)
rval_cpu, tc = run_conv_nnet2_classif(shared, isize, ksize, bsize, n_iter)
print "cpu:", rval_cpu
print "gpu:", rval_gpu
print "abs diff:", numpy.absolute(rval_gpu-rval_cpu)
......@@ -309,16 +311,21 @@ def run_test_conv_nnet2_classif(seed, isize, ksize, bsize, ignore_error=False, g
assert numpy.allclose(rval_cpu[:2], rval_gpu[:2],rtol=1e-4,atol=1e-6)
def test_lenet_28(): #MNIST
run_test_conv_nnet2_classif(23485, 28, 5, 60)
cmp_run_conv_nnet2_classif(23485, 28, 5, 60, n_iter=3)
def test_lenet_32(): #CIFAR10 / Shapeset
run_test_conv_nnet2_classif(23485, 32, 5, 60, ignore_error=False)
cmp_run_conv_nnet2_classif(23485, 32, 5, 60, ignore_error=False, n_iter=3)
def test_lenet_32_long(): #CIFAR10 / Shapeset
# this tests the gradient of downsample on the GPU,
# which does not recieve specific testing
cmp_run_conv_nnet2_classif(23485, 32, 5, 30, ignore_error=False, n_iter=50)
def test_lenet_64(): # ???
run_test_conv_nnet2_classif(23485, 64, 7, 10, ignore_error=True)
cmp_run_conv_nnet2_classif(23485, 64, 7, 10, ignore_error=False, n_iter=3)
def test_lenet_108(): # NORB
run_test_conv_nnet2_classif(23485, 108, 7, 10)
#def test_lenet_108(): # NORB
#cmp_run_conv_nnet2_classif(23485, 108, 7, 10)
def test_lenet_256(): # ImageNet
run_test_conv_nnet2_classif(23485, 256, 9, 2)
#def test_lenet_256(): # ImageNet
#cmp_run_conv_nnet2_classif(23485, 256, 9, 2)
var.py
浏览文件 @ 2413550e
......@@ -54,13 +54,16 @@ class CudaNdarraySharedVariable(SharedVariable, _operators):
if (other.type.dtype != self.dtype):
raise TypeError('Incompatible dtype', (self.dtype, other.type.dtype))
if (other.type.broadcastable != self.broadcastable):
raise TypeError('Incompatible broadcastable', (self.broadcastable, other.type.broadcastable))
raise TypeError('Incompatible broadcastable', (self, (self.broadcastable,
other.type.broadcastable)))
return GpuFromHost()(other)
CudaNdarrayType.SharedVariable = CudaNdarraySharedVariable
def shared_constructor(value, name, strict=False):
def shared_constructor(value, name, strict=False, broadcastable=None):
"""SharedVariable Constructor for TensorType"""
#TODO: what should strict mean in this context, since we always have to make a copy?
if strict:
_value = value
else:
......@@ -71,8 +74,9 @@ def shared_constructor(value, name, strict=False):
if _value.dtype.num != CudaNdarrayType.typenum:
raise TypeError('float32 ndarray required')
bcast = [0 for b in value.shape]
type = CudaNdarrayType(broadcastable=bcast)
if broadcastable is None:
broadcastable = [b==1 for b in value.shape]
type = CudaNdarrayType(broadcastable=broadcastable)
return CudaNdarraySharedVariable(type=type, value=_value, name=name, strict=strict)
......
Markdown 格式
0% 或
您添加了 0 人 到此讨论。请谨慎行事。
请先完成此评论的编辑！
请 注册 或者 后发表评论