RBFTests.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Interpolation based on Radial Basis Functions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Parameter Matrix (2D meshgrid)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2×9 Array{Int64,2}:\n",
       " 1  2  3  1  2  3  1  2  3\n",
       " 1  1  1  2  2  2  3  3  3"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P=[1 2 3 1 2 3 1 2 3;1 1 1 2 2 2 3 3 3]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "2D coordinates of the 3 points to be interpolated"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2×3 Array{Float64,2}:\n",
       " 1.5  1.1  2.5\n",
       " 1.5  1.2  2.5"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    " Pint=[1.5 1.1 2.5;1.5 1.2 2.5]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Function used to generate the response for the previous 9 datapoints "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "f (generic function with 1 method)"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f(x)=[sqrt(x[1]), x[2]^2, x[1]+x[2]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Generation of the data matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The Data - Snapshot matrix coming out from the DOE Y=\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "3×9 Array{Float64,2}:\n",
       " 1.0  1.41  1.73  1.0  1.41  1.73  1.0  1.41  1.73\n",
       " 1.0  1.0   1.0   4.0  4.0   4.0   9.0  9.0   9.0\n",
       " 2.0  3.0   4.0   3.0  4.0   5.0   4.0  5.0   6.0"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Y=zeros(3,9)\n",
    "for i in 1:9\n",
    "    Y[:,i]=f(P[:,i])\n",
    "end\n",
    "println(\"\\rThe Data - Snapshot matrix coming out from the DOE Y=\\r\")\n",
    "round.(Y;digits=2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Choice of two types of Radial Basis Function :\n",
    "## 1. Linear Spline (lin)\n",
    "## 2. Cubic Spline (cub)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. RBF : Linear Spline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "GlinMtx (generic function with 1 method)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function GlinMtx(p)\n",
    "    M = size(p, 1)\n",
    "    N = size(p, 2)\n",
    "    G = zeros(N,N)\n",
    "    for i = 1:N\n",
    "        for j = 1:N\n",
    "            G[i,j] = sqrt(sum((p[:,i] - p[:,j]).^2)) # G1[i,j]=sqrt((P[:,i][k]-P[:,j][1])^2+(P[:,i][k]-P[:,j][2])^2) for k={1:2}\n",
    "        end\n",
    "    end\n",
    "    return G\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9×9 Array{Float64,2}:\n",
       " 0.0      1.0      2.0      1.0      …  2.23607  2.0      2.23607  2.82843\n",
       " 1.0      0.0      1.0      1.41421     1.41421  2.23607  2.0      2.23607\n",
       " 2.0      1.0      0.0      2.23607     1.0      2.82843  2.23607  2.0\n",
       " 1.0      1.41421  2.23607  0.0         2.0      1.0      1.41421  2.23607\n",
       " 1.41421  1.0      1.41421  1.0         1.0      1.41421  1.0      1.41421\n",
       " 2.23607  1.41421  1.0      2.0      …  0.0      2.23607  1.41421  1.0\n",
       " 2.0      2.23607  2.82843  1.0         2.23607  0.0      1.0      2.0\n",
       " 2.23607  2.0      2.23607  1.41421     1.41421  1.0      0.0      1.0\n",
       " 2.82843  2.23607  2.0      2.23607     1.0      2.0      1.0      0.0"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Glin=GlinMtx(P)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculated Matrix A of coefficients of interpolation "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×9 Array{Float64,2}:\n",
       " 0.367804  -0.0436572  0.152806  …   0.367804  -0.0436572   0.152806\n",
       " 1.72823    0.161643   1.72823      -0.621314  -0.934163   -0.621314\n",
       " 1.30171    0.126765   0.714327      0.714327  -0.147186    0.126941"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Alin=Y*(inv(Glin'))'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "2. RBF : Cubic Spline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "GcubMtx (generic function with 1 method)"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "function GcubMtx(p)\n",
    "    M = size(p, 1)\n",
    "    N = size(p, 2)\n",
    "    G = zeros(N,N)\n",
    "    for i = 1:N\n",
    "        for j = 1:N\n",
    "            G[i,j] = sqrt(sum((p[:,i] - p[:,j]).^2))^3 # G1[i,j]=sqrt((P[:,i][k]-P[:,j][1])^2+(P[:,i][k]-P[:,j][2])^2)^3 for k={1:2}\n",
    "        end\n",
    "    end\n",
    "    return G\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9×9 Array{Float64,2}:\n",
       "  0.0       1.0       8.0       1.0      …   8.0      11.1803   22.6274\n",
       "  1.0       0.0       1.0       2.82843     11.1803    8.0      11.1803\n",
       "  8.0       1.0       0.0      11.1803      22.6274   11.1803    8.0\n",
       "  1.0       2.82843  11.1803    0.0          1.0       2.82843  11.1803\n",
       "  2.82843   1.0       2.82843   1.0          2.82843   1.0       2.82843\n",
       " 11.1803    2.82843   1.0       8.0      …  11.1803    2.82843   1.0\n",
       "  8.0      11.1803   22.6274    1.0          0.0       1.0       8.0\n",
       " 11.1803    8.0      11.1803    2.82843      1.0       0.0       1.0\n",
       " 22.6274   11.1803    8.0      11.1803       8.0       1.0       0.0"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gcub=GcubMtx(P)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Calculated Matrix A of coefficients of interpolation "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×9 Array{Float64,2}:\n",
       " 0.27631   -0.331726  0.215527  -0.383887  …  0.27631   -0.331726  0.215527\n",
       " 0.980982  -0.967142  0.980982  -1.04863      0.316732  -0.276571  0.316732\n",
       " 0.856022  -1.03782   0.68996   -1.03782      0.68996   -0.865181  0.523898"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Acub=Y*(inv(Gcub'))'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Calculate the RBF interpolated values of the points Pint to be interpolated  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. RBF : Linear Spline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "g=zeros(9,3);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Wrong\n",
    "for j in 1:3\n",
    "    for i=1:9\n",
    "    g[i,j]=sqrt((Pint[:,j][1]-P[:,j][1])^2+(Pint[:,j][2]-P[:,j][2])^2)\n",
    "end\n",
    "end\n",
    "g"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9-element Array{Float64,1}:\n",
       " 0.7071067811865476\n",
       " 0.7071067811865476\n",
       " 1.5811388300841898\n",
       " 0.7071067811865476\n",
       " 0.7071067811865476\n",
       " 1.5811388300841898\n",
       " 1.5811388300841898\n",
       " 1.5811388300841898\n",
       " 2.1213203435596424"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Node 1 (1.5,1.5)\n",
    "g₁¹=zeros(9)\n",
    "for i in 1:9\n",
    "g₁¹[i]=sqrt((1.5-P[:,i][1])^2+(1.5-P[:,i][2])^2)\n",
    "end\n",
    "g₁¹"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9-element Array{Float64,1}:\n",
       " 0.22360679774997896\n",
       " 0.9219544457292886\n",
       " 1.91049731745428\n",
       " 0.806225774829855\n",
       " 1.2041594578792296\n",
       " 2.0615528128088303\n",
       " 1.8027756377319948\n",
       " 2.0124611797498106\n",
       " 2.6172504656604803"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Node 2 (1.1,1.1)\n",
    "g₁²=zeros(9)\n",
    "for i in 1:9\n",
    "g₁²[i]=sqrt((1.1-P[:,i][1])^2+(1.2-P[:,i][2])^2)\n",
    "end\n",
    "g₁²"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9-element Array{Float64,1}:\n",
       " 2.1213203435596424\n",
       " 1.5811388300841898\n",
       " 1.5811388300841898\n",
       " 1.5811388300841898\n",
       " 0.7071067811865476\n",
       " 0.7071067811865476\n",
       " 1.5811388300841898\n",
       " 0.7071067811865476\n",
       " 0.7071067811865476"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Node 3 (2.5,2.5)\n",
    "g₁³=zeros(9)\n",
    "for i in 1:9\n",
    "g₁³[i]=sqrt((2.5-P[:,i][1])^2+(2.5-P[:,i][2])^2)\n",
    "end\n",
    "g₁³"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3-element Array{Float64,1}:\n",
       " 1.138303850496928\n",
       " 2.0769803445205515\n",
       " 2.756639112959866"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₁¹=Alin*g₁¹"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3-element Array{Float64,1}:\n",
       " 1.0177206443731015\n",
       " 1.4070346853046365\n",
       " 2.204142798417518"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₁²=Alin*g₁²"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3-element Array{Float64,1}:\n",
       " 1.5299989984668725\n",
       " 6.357504525848947\n",
       " 4.896901203624063"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₁³=Alin*g₁³"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "2. RBF : Cubic Spline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9-element Array{Float64,1}:\n",
       " 0.35355339059327384\n",
       " 0.35355339059327384\n",
       " 3.952847075210475\n",
       " 0.35355339059327384\n",
       " 0.35355339059327384\n",
       " 3.952847075210475\n",
       " 3.952847075210475\n",
       " 3.952847075210475\n",
       " 9.545941546018389"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Node 1 (1.5,1.5)\n",
    "g₂¹=zeros(9)\n",
    "for i in 1:9\n",
    "g₂¹[i]=sqrt((1.5-P[:,i][1])^2+(1.5-P[:,i][2])^2)^3\n",
    "end\n",
    "g₂¹"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9-element Array{Float64,1}:\n",
       "  0.011180339887498947\n",
       "  0.7836612788698952\n",
       "  6.973315208708123\n",
       "  0.5240467536394059\n",
       "  1.7460312139248833\n",
       "  8.761599454437528\n",
       "  5.859020822628984\n",
       "  8.150467777986734\n",
       " 17.928165689774293"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Node 2 (1.1,1.1)\n",
    "g₂²=zeros(9)\n",
    "for i in 1:9\n",
    "g₂²[i]=sqrt((1.1-P[:,i][1])^2+(1.2-P[:,i][2])^2)^3\n",
    "end\n",
    "g₂²"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "9-element Array{Float64,1}:\n",
       " 9.545941546018389\n",
       " 3.952847075210475\n",
       " 3.952847075210475\n",
       " 3.952847075210475\n",
       " 0.35355339059327384\n",
       " 0.35355339059327384\n",
       " 3.952847075210475\n",
       " 0.35355339059327384\n",
       " 0.35355339059327384"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#Node 3 (2.5,2.5)\n",
    "g₂³=zeros(9)\n",
    "for i in 1:9\n",
    "g₂³[i]=sqrt((2.5-P[:,i][1])^2+(2.5-P[:,i][2])^2)^3\n",
    "end\n",
    "g₂³"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3-element Array{Float64,1}:\n",
       " 1.2921213776313536\n",
       " 2.584422378781933\n",
       " 3.2506384352761466"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₂¹=Acub*g₂¹"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3-element Array{Float64,1}:\n",
       " 1.1367423567255974\n",
       " 1.811487244689924\n",
       " 2.5925893416474644"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₂²=Acub*g₂²"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3-element Array{Float64,1}:\n",
       " 1.6234173643382126\n",
       " 6.204892250889171\n",
       " 5.060873371329757"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₂³=Acub*g₂³"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Calculated Interpolations for the 3 points in Pint using a linear spline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×3 Array{Float64,2}:\n",
       " 1.1383   1.01772  1.53\n",
       " 2.07698  1.40703  6.3575\n",
       " 2.75664  2.20414  4.8969"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₁=[f₁¹ f₁² f₁³]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Calculated Interpolations for the 3 points in Pint using a cubic spline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×3 Array{Float64,2}:\n",
       " 1.29212  1.13674  1.62342\n",
       " 2.58442  1.81149  6.20489\n",
       " 3.25064  2.59259  5.06087"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f₂=[f₂¹ f₂² f₂³]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## To be compared to the exact solution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "3×3 Array{Float64,2}:\n",
       " 1.22474  1.04881  1.58114\n",
       " 2.25     1.44     6.25\n",
       " 3.0      2.3      5.0"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "fsol=[f([1.5,1.5]) f([1.1,1.2]) f([2.5,2.5])]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Errors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "?sum"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.5.0",
   "language": "julia",
   "name": "julia-1.5"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.5.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}