pycuda-stack-gp-2.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline \n",
    "\n",
    "import numpy as np\n",
    "import numpy\n",
    "import pandas as pd\n",
    "import re\n",
    "import random as rnd\n",
    "from sklearn.metrics import accuracy_score\n",
    "import math\n",
    "import operator\n",
    "from tqdm import tqdm_notebook as tqdm\n",
    "import random\n",
    "import warnings\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pycuda.driver as cuda\n",
    "import pycuda.autoinit\n",
    "from pycuda.compiler import SourceModule"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(891, 8)"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# read the arrays\n",
    "a = [float(x) for x in open('titanic_cpp_x_old.txt').read().split()]\n",
    "cols = int(a[0])\n",
    "bdx = array(a[1:]).reshape(-1, cols)\n",
    "bdx.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(418, 8)"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a = [float(x) for x in open('titanic_cpp_test_x_old.txt').read().split()]\n",
    "cols = int(a[0])\n",
    "bdx_test = array(a[1:]).reshape(-1, cols)\n",
    "bdx_test.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(891,)"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "bdy = array([float(x) for x in open('titanic_cpp_y.txt').read().split()])\n",
    "bdy.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(418,)"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "bdy_test = array([float(x) for x in open('titanic_cpp_test_y.txt').read().split()])\n",
    "bdy_test.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(891, 8)"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "bdx.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [],
   "source": [
    "# data organization\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_insts = 19"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_rows = 1024\n",
    "num_features = bdx.shape[1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_threads_x = 16384\n",
    "num_threads_y = num_rows"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [],
   "source": [
    "prog_len = 16\n",
    "gprogs = np.zeros((num_threads_x, prog_len), dtype=np.float32)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 2. The data is a single array of 8 features, many rows, and each thread on the Y axis takes an individual row\n",
    "gdata = np.zeros((num_rows, num_features), dtype=np.float32)\n",
    "gdatay = np.zeros((num_rows, 1), dtype=np.float32)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 3. The results - for each thread, single column of results\n",
    "gresult = np.zeros((num_rows, 1*num_threads_x), dtype=np.float32)\n",
    "gscores = np.zeros((1, 1*num_threads_x), dtype=np.float32)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [],
   "source": [
    "sorted_idx = np.zeros((1,num_threads_x), dtype=np.int32)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Allocate the buffers in GPU memory\n",
    "gprogs_gpu = cuda.mem_alloc(gprogs.nbytes)\n",
    "gprogsnew_gpu = cuda.mem_alloc(gprogs.nbytes) # to hold next generation\n",
    "sorted_idx_gpu = cuda.mem_alloc(sorted_idx.nbytes)\n",
    "\n",
    "gdata_gpu = cuda.mem_alloc(gdata.nbytes)\n",
    "gdatay_gpu = cuda.mem_alloc(gdatay.nbytes)\n",
    "gresult_gpu = cuda.mem_alloc(gresult.nbytes)\n",
    "gscores_gpu = cuda.mem_alloc(gscores.nbytes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [],
   "source": [
    "threads_per_block = 16\n",
    "blocks_x = num_threads_x // threads_per_block\n",
    "blocks_y = num_threads_y // threads_per_block"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The GP kernel\n",
    "s = f\"\"\"\n",
    "#define PROGLEN {prog_len}\n",
    "#define NUM_FEATURES {num_features}\n",
    "#define NUM_THREADS_X {num_threads_x}\n",
    "\"\"\"\n",
    "\n",
    "mod = SourceModule(s + \"\"\"\n",
    "\n",
    "#include <curand_kernel.h>\n",
    "\n",
    "#define NOP 0\n",
    "#define PUSHV 1\n",
    "#define PUSHC 2\n",
    "#define ADD 3\n",
    "#define MUL 4\n",
    "#define DIV 5\n",
    "#define NEG 6\n",
    "#define MIN 7\n",
    "#define MAX 8\n",
    "#define GREATER 9\n",
    "#define LESS 10\n",
    "#define EQUAL 11\n",
    "#define SIN 12\n",
    "#define COS 13\n",
    "#define EXP 14\n",
    "#define LOG 15\n",
    "#define SQR 16\n",
    "#define SQRT 17\n",
    "#define TANH 18\n",
    "\n",
    "\n",
    "\n",
    "const int nstates = NUM_THREADS_X;\n",
    "__device__ curandState_t* states[nstates];\n",
    "\n",
    "extern \"C\" {\n",
    "\n",
    "__global__ void initkernel(int seed)\n",
    "{\n",
    "    int tidx = threadIdx.x + blockIdx.x * blockDim.x;\n",
    "\n",
    "    if (tidx < nstates) \n",
    "    {\n",
    "        curandState_t* s = new curandState_t;\n",
    "        if (s != 0) \n",
    "        {\n",
    "            curand_init(seed, tidx, 0, s);\n",
    "        }\n",
    "\n",
    "        states[tidx] = s;\n",
    "    }\n",
    "}\n",
    "\n",
    "__device__ float randfloat(int tidx)\n",
    "{\n",
    "    curandState_t s = *states[tidx];\n",
    "    float x = curand_uniform(&s);\n",
    "    *states[tidx] = s;\n",
    "    return x;\n",
    "}\n",
    "\n",
    "__device__ unsigned int randint(int tidx, int aX, int aY)\n",
    "{\n",
    "    curandState_t s = *states[tidx];\n",
    "    unsigned int x = curand(&s);\n",
    "    *states[tidx] = s;\n",
    "    \n",
    "    if (aX == aY)\n",
    "    {\n",
    "        return aX;\n",
    "    }\n",
    "    if (aX == (aY-1))\n",
    "    {\n",
    "        // for two consecutives, pick either with equal probability\n",
    "        if (randfloat(tidx) < 0.5)\n",
    "        {\n",
    "            return aX;\n",
    "        }\n",
    "        else\n",
    "        {\n",
    "            return aY;\n",
    "        }\n",
    "    }\n",
    "    return aX + (x % (aY - aX + 1));\n",
    "}\n",
    "\n",
    "\n",
    "__global__ void init_progs(float *gprogs)\n",
    "{\n",
    "    int tidx = (blockIdx.x * blockDim.x + threadIdx.x);\n",
    "    int prog_idx =  tidx * PROGLEN;\n",
    "    float* start_idx = (gprogs + prog_idx);\n",
    "    \n",
    "    for(int ip=0; ip<(PROGLEN*2); ip += 2)\n",
    "    {\n",
    "        // random instruction\n",
    "        float inst = (float)(randint(tidx, 0, TANH));\n",
    "        float operand = 0;\n",
    "        if (inst == PUSHV)\n",
    "        {\n",
    "            operand = (float)(randint(tidx, 0, NUM_FEATURES-1));\n",
    "        }\n",
    "        if (inst == PUSHC)\n",
    "        {\n",
    "            operand = (randfloat(tidx) - 0.5) * 5;\n",
    "        }\n",
    "        \n",
    "        \n",
    "        *(start_idx + ip) = inst; // instruction\n",
    "        *(start_idx + ip + 1) = operand; // operand\n",
    "    }\n",
    "}\n",
    "\n",
    "\n",
    "// crossover\n",
    "__device__ void crossover(int tidx, float* mom, float* dad, float* baby)\n",
    "{\n",
    "    unsigned int p1 = randint(tidx, 0, PROGLEN - 1) * 2;\n",
    "\tunsigned int p2 = randint(tidx, 0, PROGLEN - 1) * 2;\n",
    "    \n",
    "\twhile (p1 == p2)\n",
    "\t{\n",
    "        p1 = randint(tidx, 0, PROGLEN - 1) * 2;\n",
    "        p2 = randint(tidx, 0, PROGLEN - 1) * 2;\n",
    "\t}\n",
    "\n",
    "\tif (p1 > p2)\n",
    "\t{\n",
    "\t\tint tmp = p1;\n",
    "\t\tp1 = p2;\n",
    "\t\tp2 = tmp;\n",
    "\t}\n",
    "    \n",
    "    float x;\n",
    "    \n",
    "    for (unsigned int i = 0; i < p1; i += 2)\n",
    "\t{\n",
    "        x = mom[i];\n",
    "        baby[i] = x;\n",
    "        x =  mom[i+1];\n",
    "        baby[i+1] = x;\n",
    "\t}\n",
    "    \n",
    "\tfor (unsigned int i = p1; i < p2; i += 2)\n",
    "\t{\n",
    "        x = dad[i];\n",
    "        baby[i] = x;\n",
    "        x =  dad[i+1];\n",
    "        baby[i+1] = x;\n",
    "\t}\n",
    "    \n",
    "\tfor (unsigned int i = p2; i < PROGLEN*2; i += 2)\n",
    "\t{\n",
    "        x = mom[i];\n",
    "        baby[i] = x;\n",
    "        x =  mom[i+1];\n",
    "        baby[i+1] = x;\n",
    "\t}\n",
    "}\n",
    "\n",
    "\n",
    "// mutation\n",
    "__device__ void mutate(int tidx, float* baby)\n",
    "{\n",
    "    for (int i = 0; i < PROGLEN*2; i += 2)\n",
    "\t{\n",
    "\t\tif (randfloat(tidx) < 0.1)\n",
    "\t\t{\n",
    "\t\t\tif ((baby[i] == PUSHV) && (randfloat(tidx) < 0.8))\n",
    "\t\t\t{\n",
    "\t\t\t\tbaby[i + 1] = (float)(randint(tidx, 0, NUM_FEATURES - 1));\n",
    "\t\t\t}\n",
    "\t\t\telse if ((baby[i] == PUSHC) && (randfloat(tidx) < 0.8))\n",
    "\t\t\t{\n",
    "\t\t\t\tbaby[i + 1] = (randfloat(tidx) - 0.5) * 5;\n",
    "\t\t\t}\n",
    "\t\t\telse\n",
    "\t\t\t{\n",
    "                float inst = (float)(randint(tidx, 0, TANH));\n",
    "                float operand = 0;\n",
    "                if (inst == PUSHV)\n",
    "                {\n",
    "                    operand = (float)(randint(tidx, 0, NUM_FEATURES-1));\n",
    "                }\n",
    "                if (inst == PUSHC)\n",
    "                {\n",
    "                    operand = (randfloat(tidx) - 0.5) * 5;\n",
    "                }\n",
    "                baby[i] = inst;\n",
    "                baby[i+1] = operand;\n",
    "\t\t\t}\n",
    "\t\t}\n",
    "\t}\n",
    "}\n",
    "\n",
    "\n",
    "\n",
    "// will produce one baby per thread and store it in gprogsnew\n",
    "__global__ void reproduce(int* sorted_idx, float *gprogs, float *gprogsnew)\n",
    "{\n",
    "    int tidx = (blockIdx.x * blockDim.x + threadIdx.x);\n",
    "    int prog_idx =  tidx * PROGLEN;\n",
    "    float* baby_start_idx = (gprogsnew + prog_idx);\n",
    "    \n",
    "    // decide whether to do crossover or not\n",
    "    if (randfloat(tidx) < 0.8)\n",
    "    {\n",
    "        // will do crossover\n",
    "        int idx1 = sorted_idx[randint(tidx, 0, (int)((float)(NUM_THREADS_X) * 0.2))];\n",
    "        int idx2 = sorted_idx[randint(tidx, 0, (int)((float)(NUM_THREADS_X) * 0.2))];\n",
    "        \n",
    "        float* mom = gprogs + idx1*PROGLEN;\n",
    "        float* dad = gprogs + idx2*PROGLEN;\n",
    "        \n",
    "        crossover(tidx, mom, dad, baby_start_idx);\n",
    "        mutate(tidx, baby_start_idx);\n",
    "    }\n",
    "    else \n",
    "    {\n",
    "        // will only mutate, but copy the baby first\n",
    "        int idx = sorted_idx[randint(tidx, 0, (int)((float)(NUM_THREADS_X) * 0.2))];\n",
    "        float* mom = gprogs + idx*PROGLEN;\n",
    "        \n",
    "        for(int i=0; i<PROGLEN*2; i++)\n",
    "        {\n",
    "            baby_start_idx[i] = mom[i];\n",
    "        }\n",
    "        \n",
    "        if (randfloat(tidx) < 0.5)\n",
    "        {\n",
    "            mutate(tidx, baby_start_idx);\n",
    "        }\n",
    "    }\n",
    "}\n",
    "\n",
    "\n",
    "__global__ void compute(int num_real_rows, float *gprogs, float *gdata, float *gresult)\n",
    "{\n",
    "    int prog_idx = (blockIdx.x * blockDim.x + threadIdx.x) * PROGLEN;\n",
    "    int data_row_idx = (blockIdx.y * blockDim.y + threadIdx.y) * NUM_FEATURES;\n",
    "\n",
    "    int th_y = blockIdx.y * blockDim.y + threadIdx.y;\n",
    "    int th_x = blockIdx.x * blockDim.x + threadIdx.x;\n",
    "\n",
    "    int result_idx = NUM_THREADS_X * th_y + th_x; \n",
    "\n",
    "    if (th_y >= num_real_rows)\n",
    "        return; // don't compute where there's no data\n",
    "\n",
    "    float stack[PROGLEN+4]; // a little bigger just in case\n",
    "\n",
    "    int sp = 0;\n",
    "\n",
    "    // push a 0 to the stack\n",
    "    stack[sp] = 0;\n",
    "    sp++;\n",
    "\n",
    "    float x = 0;\n",
    "    float y = 0;\n",
    "\n",
    "    float inst = 0;\n",
    "    float operand = 0;\n",
    "\n",
    "    /////////////////////////////////\n",
    "    // Evaluate program on data row\n",
    "    float* start_idx = (gprogs + prog_idx);\n",
    "\n",
    "    for(int ip=0; ip<(PROGLEN*2); ip += 2)\n",
    "    {\n",
    "        // fetch instruction and operand\n",
    "        inst = *(start_idx + ip);\n",
    "        operand = *(start_idx + ip + 1);\n",
    "\n",
    "        // execute instruction\n",
    "        if (inst == NOP)\n",
    "            continue; \n",
    "\n",
    "        if (inst == PUSHV)\n",
    "        {\n",
    "            // fetch variable from data\n",
    "            x = *(gdata + data_row_idx + (int)(operand));\n",
    "\n",
    "            // push to stack\n",
    "            stack[sp] = x;\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if (inst == PUSHC)\n",
    "        {\n",
    "            // push constant to stack\n",
    "            stack[sp] = operand;\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        // math/logic instructions go below\n",
    "        if ((inst == ADD) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = x + y;\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == MUL) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = x * y;\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == DIV) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = x / y;\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == NEG) && (sp >= 1))\n",
    "        {\n",
    "            /*\n",
    "            // pop one value from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = -x;\n",
    "            sp++;*/\n",
    "\n",
    "            stack[sp] = -stack[sp];\n",
    "        }\n",
    "\n",
    "        if ((inst == MIN) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = (x > y)?y:x;\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == MAX) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = (x > y)?x:y;\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == GREATER) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = (float)(x > y);\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == LESS) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = (float)(x < y);\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == EQUAL) && (sp >= 2))\n",
    "        {\n",
    "            // pop two values from stack\n",
    "            x = stack[sp];\n",
    "            sp--;\n",
    "            y = stack[sp];\n",
    "            sp--;\n",
    "\n",
    "            // push result to stack\n",
    "            stack[sp] = (float)(x == y);\n",
    "            sp++;\n",
    "        }\n",
    "\n",
    "        if ((inst == SIN) && (sp >= 1))\n",
    "        {\n",
    "            stack[sp] = sin(stack[sp]);\n",
    "        }\n",
    "\n",
    "        if ((inst == COS) && (sp >= 1))\n",
    "        {\n",
    "            stack[sp] = cos(stack[sp]);\n",
    "        }\n",
    "\n",
    "        if ((inst == EXP) && (sp >= 1))\n",
    "        {\n",
    "            stack[sp] = exp(stack[sp]);\n",
    "        }\n",
    "\n",
    "        if ((inst == LOG) && (sp >= 1))\n",
    "        {\n",
    "            stack[sp] = log(stack[sp]);\n",
    "        }\n",
    "\n",
    "        if ((inst == SQR) && (sp >= 1))\n",
    "        {\n",
    "            stack[sp] = stack[sp]*stack[sp];\n",
    "        }\n",
    "\n",
    "        if ((inst == SQRT) && (sp >= 1))\n",
    "        {\n",
    "            stack[sp] = sqrt(stack[sp]);\n",
    "        }\n",
    "\n",
    "        if ((inst == TANH) && (sp >= 1))\n",
    "        {\n",
    "            stack[sp] = tanh(stack[sp]);\n",
    "        }\n",
    "    }\n",
    "\n",
    "    // store result\n",
    "    float f = 1/(1+exp(-stack[sp]));\n",
    "    if (f > 0.5) \n",
    "    {\n",
    "        f = 1.0;\n",
    "    }\n",
    "    else \n",
    "    {\n",
    "        f = 0.0;\n",
    "    }\n",
    "    gresult[result_idx] = f;    \n",
    "}\n",
    "\n",
    "}\n",
    "  \"\"\", no_extern_c=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [],
   "source": [
    "func = mod.get_function(\"compute\")\n",
    "init_func = mod.get_function(\"initkernel\")\n",
    "initprog_func = mod.get_function(\"init_progs\")\n",
    "reproduce_func = mod.get_function(\"reproduce\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 48,
   "metadata": {},
   "outputs": [],
   "source": [
    "modsc = SourceModule(\"\"\"\n",
    "\n",
    "#define NUM_THREADS_X %d\n",
    "\n",
    "  __global__ void scores(int num_real_rows, float *gdatay, float *gresult, float *gscores)\n",
    "  {\n",
    "      int th_x = blockIdx.x * blockDim.x + threadIdx.x;\n",
    "      float sum = 0;\n",
    "      \n",
    "      for(int i=0; i<num_real_rows; i++)\n",
    "      {\n",
    "          float y1 = gdatay[i];\n",
    "          float y2 = *(gresult + i*NUM_THREADS_X + th_x);\n",
    "          if (y1 == y2)\n",
    "          {\n",
    "              sum ++;\n",
    "          }\n",
    "      }\n",
    "      \n",
    "      sum /= num_real_rows;\n",
    "      gscores[th_x] = sum;\n",
    "  }\n",
    "  \"\"\" % (num_threads_x))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [],
   "source": [
    "funcsc = modsc.get_function(\"scores\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 50,
   "metadata": {},
   "outputs": [],
   "source": [
    "# call this only once\n",
    "init_func(np.int32(123456789), block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 51,
   "metadata": {},
   "outputs": [],
   "source": [
    "initprog_func(gprogs_gpu, block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 105,
   "metadata": {},
   "outputs": [],
   "source": [
    "cuda.memcpy_dtoh(gprogs, gprogs_gpu)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 106,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[10.,  0.,  5., ...,  0.,  0.,  0.],\n",
       "       [15.,  0., 13., ...,  0., 13.,  0.],\n",
       "       [10.,  0.,  3., ...,  0., 13.,  0.],\n",
       "       ...,\n",
       "       [ 9.,  0.,  1., ...,  0., 13.,  0.],\n",
       "       [ 6.,  0.,  5., ...,  0., 11.,  0.],\n",
       "       [10.,  0., 10., ...,  0., 13.,  0.]], dtype=float32)"
      ]
     },
     "execution_count": 106,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gprogs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 100,
   "metadata": {},
   "outputs": [],
   "source": [
    "# one evaluattion to get the scores\n",
    "gdata[0 : bdx.shape[0], :] = bdx\n",
    "gdatay[0 : bdx.shape[0], :] = bdy.reshape(-1, 1)\n",
    "\n",
    "# random search\n",
    "num_real_rows = np.int32(bdx.shape[0])\n",
    "\n",
    "cuda.memcpy_htod(gdata_gpu, gdata)\n",
    "cuda.memcpy_htod(gdatay_gpu, gdatay)\n",
    "\n",
    "func(num_real_rows, gprogs_gpu, gdata_gpu, gresult_gpu, \n",
    "         block=(threads_per_block, threads_per_block, 1), grid=(blocks_x, blocks_y, 1))\n",
    "    \n",
    "funcsc(num_real_rows, gdatay_gpu, gresult_gpu, gscores_gpu, \n",
    " block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))\n",
    "\n",
    "cuda.memcpy_dtoh(gscores, gscores_gpu)\n",
    "accs = gscores.reshape(-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 101,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([0.67901236, 0.61616164, 0.3838384 , ..., 0.61616164, 0.6386083 ,\n",
       "       0.6386083 ], dtype=float32)"
      ]
     },
     "execution_count": 101,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "accs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 102,
   "metadata": {},
   "outputs": [],
   "source": [
    "sorted_idx = np.argsort(accs)[::-1].astype(np.int32)\n",
    "# copy that to the GPU\n",
    "cuda.memcpy_htod(sorted_idx_gpu, sorted_idx)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 103,
   "metadata": {},
   "outputs": [],
   "source": [
    "reproduce_func(sorted_idx_gpu, gprogs_gpu, gprogsnew_gpu, block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 104,
   "metadata": {},
   "outputs": [],
   "source": [
    "# copy the new pop over the old pop\n",
    "cuda.memcpy_dtod(gprogs_gpu, gprogsnew_gpu, gprogs.nbytes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gdata[0 : bdx.shape[0], :] = bdx\n",
    "gdatay[0 : bdx.shape[0], :] = bdy.reshape(-1, 1)\n",
    "\n",
    "cuda.memcpy_htod(gdata_gpu, gdata)\n",
    "cuda.memcpy_htod(gdatay_gpu, gdatay)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "best_ever = 0\n",
    "for i in tqdm(range(10000000)):\n",
    "    \n",
    "    #cuda.memcpy_htod(gprogs_gpu, gprogs)\n",
    "    initprog_func(gprogs_gpu, block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))\n",
    "\n",
    "    func(num_real_rows, gprogs_gpu, gdata_gpu, gresult_gpu, \n",
    "         block=(threads_per_block, threads_per_block, 1), grid=(blocks_x, blocks_y, 1))\n",
    "    \n",
    "    funcsc(num_real_rows, gdatay_gpu, gresult_gpu, gscores_gpu, \n",
    "     block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))\n",
    "\n",
    "    cuda.memcpy_dtoh(gscores, gscores_gpu)\n",
    "    accs = gscores.reshape(-1)\n",
    "    best = np.max(accs)\n",
    "    if best > best_ever:\n",
    "        best_ever = best\n",
    "        print('new best', best)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this will hold the result\n",
    "gresult = np.zeros((num_rows, 1*num_threads_x), dtype=np.float32)\n",
    "gresult, gresult.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# this will hold the scores\n",
    "gscores = np.zeros((1, 1*num_threads_x), dtype=np.float32)\n",
    "gscores, gscores.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# copy arrays to GPU\n",
    "cuda.memcpy_htod(gprogs_gpu, gprogs)\n",
    "\n",
    "cuda.memcpy_htod(gdata_gpu, gdata)\n",
    "cuda.memcpy_htod(gdatay_gpu, gdatay)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "num_real_rows = np.array(bdx.shape).astype(np.int32)[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# execute kernel\n",
    "func(num_real_rows, gprogs_gpu, gdata_gpu, gresult_gpu, \n",
    "     block=(threads_per_block, threads_per_block, 1), grid=(blocks_x, blocks_y, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# compute the scores now\n",
    "funcsc(num_real_rows, gdatay_gpu, gresult_gpu, gscores_gpu, \n",
    "     block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# get back the data and display it\n",
    "cuda.memcpy_dtoh(gscores, gscores_gpu)\n",
    "gscores.reshape(-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# sort the arrays by score\n",
    "hs_idx = np.argsort(gscores.reshape(-1))[::-1]\n",
    "gprogs = gprogs[hs_idx]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# choose parent indices\n",
    "pr_idx = np.random.randint(0, int(len(hs_idx)*0.2), size=(len(hs_idx),))\n",
    "pr_idx"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# replace pop with the new individuals, but flatten them\n",
    "gprogs = gprogs[pr_idx].reshape(-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# mutate the pop randomly\n",
    "numtomut = int( len(pr_idx)*0.3 )\n",
    "mt_idx = np.random.randint(0, numtomut, size=(numtomut,))\n",
    "\n",
    "gprogs[mt_idx] = np.random.randint(0, num_insts, size=(numtomut,)).astype(dtype=np.float32)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gprogs = gprogs.reshape(len(hs_idx),-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gprogs = np.zeros((num_threads_x, prog_len), dtype=np.float32)\n",
    "inst_idx = np.arange(0, prog_len, 2)\n",
    "gprogs[:,inst_idx] = np.random.randint(0, num_insts, size=(num_threads_x, prog_len//2))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# %%timeit \n",
    "\n",
    "# gresult = np.zeros((num_rows, 1*num_threads_x), dtype=np.float32)\n",
    "\n",
    "# cuda.memcpy_htod(gops_gpu, gops)\n",
    "# cuda.memcpy_htod(gvars_gpu, gvars)\n",
    "# cuda.memcpy_htod(gconsts_gpu, gconsts)\n",
    "# cuda.memcpy_htod(gis_var_gpu, gis_var)\n",
    "# cuda.memcpy_htod(gdata_gpu, gdata)\n",
    "\n",
    "# func(num_real_rows, gops_gpu, gvars_gpu, gconsts_gpu, gis_var_gpu, gdata_gpu, gresult_gpu, \n",
    "#      block=(threads_per_block, threads_per_block, 1), grid=(blocks_x, blocks_y, 1))\n",
    "\n",
    "# funcsc(num_real_rows, gdatay_gpu, gresult_gpu, gscores_gpu, \n",
    "#      block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))\n",
    "\n",
    "# cuda.memcpy_dtoh(gscores, gscores_gpu)\n",
    "\n",
    "# # new generation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "hist(gscores.reshape(-1), 50);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cuda.memcpy_htod(gdata_gpu, gdata)\n",
    "cuda.memcpy_htod(gdatay_gpu, gdatay)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# initial random population \n",
    "\n",
    "# gprogs = np.zeros((num_threads_x, prog_len), dtype=np.float32)\n",
    "# inst_idx = np.arange(0, prog_len, 2)\n",
    "# gprogs[:,inst_idx] = np.random.randint(0, num_insts, size=(num_threads_x, prog_len//2))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "best_ever = 0\n",
    "for i in tqdm(range(10000000)):\n",
    "    \n",
    "    #cuda.memcpy_htod(gprogs_gpu, gprogs)\n",
    "\n",
    "    func(num_real_rows, gprogs_gpu, gdata_gpu, gresult_gpu, \n",
    "         block=(threads_per_block, threads_per_block, 1), grid=(blocks_x, blocks_y, 1))\n",
    "    \n",
    "    funcsc(num_real_rows, gdatay_gpu, gresult_gpu, gscores_gpu, \n",
    "     block=(threads_per_block, 1, 1), grid=(blocks_x, 1, 1))\n",
    "\n",
    "    cuda.memcpy_dtoh(gscores, gscores_gpu)\n",
    "    accs = gscores.reshape(-1)\n",
    "    best = np.max(accs)\n",
    "    if best > best_ever:\n",
    "        best_ever = best\n",
    "        print('new best', best)\n",
    "        \n",
    "   "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.8"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}