Metavoxel GPU Validation 2.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "using LinearAlgebra\n",
    "using Plots\n",
    "import JSON\n",
    "# using Quaternions\n",
    "using StaticArrays, Rotations\n",
    "using Distributed\n",
    "using StaticArrays, BenchmarkTools\n",
    "using Base.Threads\n",
    "using CUDAnative\n",
    "using CuArrays,CUDAdrv \n",
    "using Test\n",
    "import Base: +, * , -, ^"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Todo\n",
    "- create struct for material and get for each edge its properties\n",
    "- implement getTimestep\n",
    "- implement on single voxel\n",
    "- get reat E and L\n",
    "- compare to Frame3dd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 54,
   "metadata": {},
   "outputs": [],
   "source": [
    "struct Vector3\n",
    "    x::Float64\n",
    "    y::Float64\n",
    "    z::Float64\n",
    "    function Vector3()\n",
    "        x=0.0\n",
    "        y=0.0\n",
    "        z=0.0\n",
    "        new(x,y,z)\n",
    "    end\n",
    "    function Vector3(x,y,z)\n",
    "       new(x,y,z)\n",
    "    end\n",
    "end\n",
    "struct Quaternion\n",
    "    x::Float64\n",
    "    y::Float64\n",
    "    z::Float64\n",
    "    w::Float64\n",
    "    function Quaternion()\n",
    "        x=0.0\n",
    "        y=0.0\n",
    "        z=0.0\n",
    "        w=1.0\n",
    "        new(x,y,z,w)\n",
    "    end\n",
    "    function Quaternion(x,y,z,w)\n",
    "        new(x,y,z,w)\n",
    "    end\n",
    "end\n",
    "struct RotationMatrix\n",
    "    te1::Float64\n",
    "    te2::Float64\n",
    "    te3::Float64\n",
    "    te4::Float64\n",
    "    te5::Float64\n",
    "    te6::Float64\n",
    "    te7::Float64\n",
    "    te8::Float64\n",
    "    te9::Float64\n",
    "    te10::Float64\n",
    "    te11::Float64\n",
    "    te12::Float64\n",
    "    te13::Float64\n",
    "    te14::Float64\n",
    "    te15::Float64\n",
    "    te16::Float64\n",
    "    function RotationMatrix()\n",
    "        te1 =0.0\n",
    "        te2 =0.0\n",
    "        te3 =0.0\n",
    "        te4 =0.0\n",
    "        te5 =0.0\n",
    "        te6 =0.0\n",
    "        te7 =0.0\n",
    "        te8 =0.0\n",
    "        te9 =0.0\n",
    "        te10=0.0\n",
    "        te11=0.0\n",
    "        te12=0.0\n",
    "        te13=0.0\n",
    "        te14=0.0\n",
    "        te15=0.0\n",
    "        te16=0.0\n",
    "        new(te1,te2,te3,te4,te5,te6,te7,te8,te9,te10,te11,te12,te13,te14,te15,te16)\n",
    "    end\n",
    "    function RotationMatrix(te1,te2,te3,te4,te5,te6,te7,te8,te9,te10,te11,te12,te13,te14,te15,te16)\n",
    "        new(te1,te2,te3,te4,te5,te6,te7,te8,te9,te10,te11,te12,te13,te14,te15,te16)\n",
    "    end\n",
    "end\n",
    "\n",
    "+(f::Vector3, g::Vector3)=Vector3(f.x+g.x , f.y+g.y,f.z+g.z )\n",
    "-(f::Vector3, g::Vector3)=Vector3(f.x-g.x , f.y-g.y,f.z-g.z )\n",
    "*(f::Vector3, g::Vector3)=Vector3(f.x*g.x , f.y*g.y,f.z*g.z )\n",
    "\n",
    "+(f::Vector3, g::Number)=Vector3(f.x+g , f.y+g,f.z+g )\n",
    "-(f::Vector3, g::Number)=Vector3(f.x-g , f.y-g,f.z-g )\n",
    "*(f::Vector3, g::Number)=Vector3(f.x*g , f.y*g,f.z*g )\n",
    "\n",
    "+(g::Vector3, f::Number)=Vector3(f.x+g , f.y+g,f.z+g )\n",
    "-(g::Vector3, f::Number)=Vector3(g-f.x , g-f.y,g-f.z )\n",
    "*(g::Vector3, f::Number)=Vector3(f.x*g , f.y*g,f.z*g )\n",
    "\n",
    "addX(f::Vector3, g::Number)=Vector3(f.x+g , f.y,f.z)\n",
    "addY(f::Vector3, g::Number)=Vector3(f.x , f.y+g,f.z )\n",
    "addZ(f::Vector3, g::Number)=Vector3(f.x , f.y,f.z+g )\n",
    "\n",
    "function normalizeVector3(f::Vector3)\n",
    "    leng=sqrt((f.x * f.x) + (f.y * f.y) + (f.z * f.z))\n",
    "    return Vector3(f.x/leng,f.y/leng,f.z/leng)\n",
    "    \n",
    "end\n",
    "function normalizeQuaternion(f::Quaternion)\n",
    "    l = sqrt((f.x * f.x) + (f.y * f.y) + (f.z * f.z)+ (f.w * f.w))\n",
    "    if l === 0 \n",
    "        qx = 0\n",
    "        qy = 0\n",
    "        qz = 0\n",
    "        qw = 1\n",
    "    else \n",
    "        l = 1 / l\n",
    "        qx = f.x * l\n",
    "        qy = f.y * l\n",
    "        qz = f.z * l\n",
    "        qw = f.w * l\n",
    "    end\n",
    "    return Quaternion(qx,qy,qz,qw)\n",
    "end\n",
    "\n",
    "function normalizeQuaternion1!(fx::Float64,fy::Float64,fz::Float64,fw::Float64)\n",
    "    l = sqrt((fx * fx) + (fy * fy) + (fz * fz)+ (fw * fw))\n",
    "    if l === 0 \n",
    "        qx = 0.0\n",
    "        qy = 0.0\n",
    "        qz = 0.0\n",
    "        qw = 1.0\n",
    "    else \n",
    "        l = 1.0 / l\n",
    "        qx = fx * l\n",
    "        qy = fy * l\n",
    "        qz = fz * l\n",
    "        qw = fw * l\n",
    "    end\n",
    "    return qx,qy,qz,qw\n",
    "end\n",
    "\n",
    "\n",
    "function dotVector3(f::Vector3, g::Vector3)\n",
    "    return (f.x * g.x) + (f.y * g.y) + (f.z * g.z)\n",
    "end\n",
    "\n",
    "function Base.show(io::IO, v::Vector3)\n",
    "    print(io, \"x:$(v.x), y:$(v.y), z:$(v.z)\")\n",
    "end\n",
    "\n",
    "function Base.show(io::IO, v::Quaternion)\n",
    "    print(io, \"x:$(v.x), y:$(v.y), z:$(v.z), w:$(v.z)\")\n",
    "end\n",
    "\n",
    "Base.Broadcast.broadcastable(q::Vector3) = Ref(q)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 55,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "simulateParallel (generic function with 2 methods)"
      ]
     },
     "execution_count": 55,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function simulateParallel(numTimeSteps,dt)\n",