SHARE
TWEET

Untitled

a guest Jul 8th, 2019 6 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. import numpy
  2. import theano
  3. from theano import gof
  4. import theano.tensor as T
  5.  
  6. theano.config.gcc.cxxflags="<NEED TO ENTER PATH TO INCLUDES HERE>"
  7. class wall_log_op(gof.Op):
  8.     __props__ = ()
  9.  
  10.     def make_node(self, value, n, m, N, odds):
  11.         value_ =T.as_tensor_variable(numpy.array(value, dtype=numpy.int32))
  12.         n_ =T.as_tensor_variable(numpy.array(n, dtype=numpy.int32))
  13.         m_ =T.as_tensor_variable(numpy.array(m, dtype=numpy.int32))
  14.         N_ =T.as_tensor_variable(numpy.array(N, dtype=numpy.int32))
  15.         odds_ = T.as_tensor_variable(numpy.array(odds, dtype=numpy.double))
  16.  
  17.         output_var = odds_.type()
  18.  
  19.         return gof.Apply(self, [value_, n_, m_, N_, odds_], [output_var])
  20.  
  21.     def c_code_cache_version(self):
  22.         return (1, 0)
  23.  
  24.     def c_headers(self):
  25.         return (
  26.         ["<time.h>", "<randomc.h>", "<mersenne.cpp>", "<stocc.h>", "<stoc1.cpp>", "<stoc3.cpp>", "<userintf.cpp>",
  27.          "<wnchyppr.cpp>", "<fnchyppr.cpp>"])
  28.  
  29.     def c_code(self, node, name, inp, out, sub):
  30.  
  31.         value, n, m, N, odds = inp
  32.         logp, = out
  33.  
  34.         # Extract the dtypes of the inputs and outputs storage to
  35.         # be able to declare pointers for those dtypes in the C
  36.         # code.
  37.         dtype_value = node.inputs[0].dtype
  38.  
  39.         dtype_n = node.inputs[1].dtype
  40.         dtype_m = node.inputs[2].dtype
  41.         dtype_N = node.inputs[3].dtype
  42.         dtype_odds = node.inputs[4].dtype
  43.         dtype_logp = node.outputs[0].dtype
  44.  
  45.        
  46.         itemsize_value = numpy.dtype(dtype_value).itemsize
  47.    
  48.         itemsize_n = numpy.dtype(dtype_n).itemsize
  49.         itemsize_m = numpy.dtype(dtype_m).itemsize
  50.         itemsize_N = numpy.dtype(dtype_N).itemsize
  51.         itemsize_odds = numpy.dtype(dtype_odds).itemsize
  52.         itemsize_logp = numpy.dtype(dtype_logp).itemsize
  53.  
  54.         fail = sub['fail']
  55.  
  56.         c_code = """
  57.  
  58.        // Validate that the output storage exists and has the same
  59.        // dimension as x.
  60.        if (NULL == %(logp)s  ||
  61.            PyArray_DIMS(%(logp)s)[0] != PyArray_DIMS(%(value)s)[0])
  62.        {
  63.  
  64.            /* Reference received to invalid output variable.
  65.            Decrease received reference's ref count and allocate new
  66.            output variable */
  67.  
  68.            Py_XDECREF(%(logp)s);
  69.            %(logp)s = (PyArrayObject*)PyArray_EMPTY(1,
  70.                                                PyArray_DIMS(%(value)s),
  71.                                                PyArray_TYPE(%(odds)s),
  72.                                                0);
  73.            if (!%(logp)s) {
  74.                %(fail)s;
  75.            }
  76.      
  77.        }
  78.  
  79.  
  80.        {
  81.  
  82.            npy_%(dtype_logp)s* logp_data_ptr =
  83.                            (npy_%(dtype_logp)s*)PyArray_DATA(%(logp)s);
  84.  
  85.            npy_%(dtype_value)s* value_data_ptr =
  86.                            (npy_%(dtype_value)s*)PyArray_DATA(%(value)s);
  87.            npy_%(dtype_n)s* n_data_ptr =
  88.                            (npy_%(dtype_n)s*)PyArray_DATA(%(n)s);
  89.  
  90.            npy_%(dtype_m)s* m_data_ptr =
  91.                            (npy_%(dtype_m)s*)PyArray_DATA(%(m)s);
  92.  
  93.            npy_%(dtype_N)s* N_data_ptr =
  94.                            (npy_%(dtype_N)s*)PyArray_DATA(%(N)s);
  95.  
  96.            npy_%(dtype_odds)s* odds_data_ptr =
  97.                            (npy_%(dtype_odds)s*)PyArray_DATA(%(odds)s);
  98.  
  99.  
  100.            int value_stride = PyArray_STRIDES(%(value)s)[0] / %(itemsize_value)s;
  101.  
  102.            int n_stride = PyArray_STRIDES(%(n)s)[0] / %(itemsize_n)s;
  103.  
  104.            int m_stride = PyArray_STRIDES(%(m)s)[0] / %(itemsize_m)s;
  105.            int N_stride = PyArray_STRIDES(%(N)s)[0] / %(itemsize_N)s;
  106.            int odds_stride = PyArray_STRIDES(%(odds)s)[0] / %(itemsize_odds)s;
  107.            int logp_stride = PyArray_STRIDES(%(logp)s)[0] / %(itemsize_logp)s;
  108.            int value_dim = PyArray_DIMS(%(value)s)[0];
  109.  
  110.  
  111.  
  112.            for(int i=0; i < value_dim; i++)
  113.            {
  114.             CWalleniusNCHypergeometric wnch(
  115.                 n_data_ptr[i * n_stride],
  116.                 m_data_ptr[i * m_stride],
  117.                 N_data_ptr[i * N_stride],
  118.                 odds_data_ptr[i * odds_stride],
  119.                                1.E-8);
  120.                logp_data_ptr[i * logp_stride] =  (log(wnch.probability(value_data_ptr[i * value_stride])));
  121.  
  122.            }
  123.  
  124.            
  125.        }
  126.        """
  127.  
  128.         return c_code % locals()
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!
 
Top