Closed_Form_CF.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Final code"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## [Link for CF Solution](https://www.overleaf.com/4679963215kkhxcqqhzrxh)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Analysis based on delta coefficients"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from pylab import rcParams\n",
    "import pandas as pd\n",
    "import os"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Parametres"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "N2, NR, N1, P1 = 1.0, 1.0, 10.0, 10.0\n",
    "P2_max, PR_max = 10.0, 10.0\n",
    "tau = 0.25"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Functions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_A(H11): \n",
    "    return ((((H11**2)*P1)/(((1+(((H11**2)*P1)/(N1)))**(1-tau))-1))-N1)    \n",
    "\n",
    "def C(x): \n",
    "    return (1/2*np.log2(1+x))\n",
    "\n",
    "def mult_res_v2(F, hR2, h12, h1R, hR1, h22, h2R, h21, h11):\n",
    "  # calculer BF ou SA pour plusieurs echantillons \n",
    "  all_res = np.zeros((hR2.shape[0], 3))\n",
    "  for i in range (hR2.shape[0]):\n",
    "        \n",
    "        all_res[i][0] += F(hR2[i], h12[i], h1R[i], hR1[i], h22[i], h2R[i], h21[i], h11[i])[0]\n",
    "        all_res[i][1] += F(hR2[i], h12[i], h1R[i], hR1[i], h22[i], h2R[i], h21[i], h11[i])[1]\n",
    "        all_res[i][2] += F(hR2[i], h12[i], h1R[i], hR1[i], h22[i], h2R[i], h21[i], h11[i])[2]\n",
    "\n",
    "  return all_res\n",
    " \n",
    "\n",
    "def Delta(HR2, H12, H1R, HR1, H22, H2R, H21, H11):\n",
    "   \n",
    "    #A = (((H11**2*P1)/(((1+((H11**2*P1)/(N1)))**(1-tau))-1))-N1)    \n",
    "    A = calculate_A(H11)\n",
    "    N2_tilde = ((H12**2)*P1+N2)\n",
    "    NR_tilde = ((H1R**2)*P1+NR)\n",
    "    P_Z = ((H12*H1R*P1)/(np.sqrt(N2_tilde*NR_tilde)))\n",
    "    K1 = ((H2R**2)*N2_tilde+(H22**2)*NR_tilde-2*H2R*H22*P_Z*np.sqrt(N2_tilde*NR_tilde))\n",
    "    K2 = ((1-P_Z**2)*NR_tilde*N2_tilde)\n",
    "    \n",
    "    C1 = K1*(HR1**2)*((H22**2)*(HR1**2)-(HR2**2)*(H21**2))\n",
    "    C2 = K1*(HR2**2)*(H21**2)*A-2*(HR1**2)*(H22**2)*K1*A-(HR1**2)*(H22**2)*(H21**2)*K2\n",
    "    C3 = (H22**2)*K1*(A**2)+(H22**2)*(H21**2)*K2*A\n",
    "    C4 = K2*(HR2**2)*(H21**2)**2-N2_tilde*K1*(H21**2)*(HR1**2)\n",
    "    C5 = N2_tilde*(H21**2)*(K1*A+K2*(H21**2))\n",
    "    D = C1**2*C5**2-C1*C4*(C2*C5-C3*C4)\n",
    "    \n",
    "    return C1,C2,C3,C4,C5,D\n",
    "\n",
    "def f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr, P2):\n",
    "    # f en fonction de PR et P2\n",
    "    A = calculate_A(H11)     \n",
    "\n",
    "    N2_tilde = ((H12**2)*P1+N2)\n",
    "    NR_tilde = ((H1R**2)*P1+NR)\n",
    "    P_Z = ((H12*H1R*P1)/(np.sqrt(N2_tilde*NR_tilde)))\n",
    "    K1 = ((H2R**2)*N2_tilde+(H22**2)*NR_tilde-2*H2R*H22*P_Z*np.sqrt(N2_tilde*NR_tilde))\n",
    "    K2 = ((1-P_Z**2)*NR_tilde*N2_tilde)\n",
    "    #res = (K1*HR2**2*P2*Pr+H22**2*P2*K1*P2+H22**2*P2*K2)/(K2*HR2**2*Pr+N2_tilde*K1*P2+N2_tilde*K2)\n",
    "    res = (K1*(HR2**2)*P2*Pr+(H22**2)*P2*(K1*P2+K2))/(K2*(HR2**2)*Pr+N2_tilde*K1*P2+N2_tilde*K2)\n",
    "\n",
    "    \n",
    "    return res\n",
    "\n",
    "def f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x):\n",
    "    'fonction depend seulment de Pr (x apres changement de variable) '\n",
    "    # f en fonction de PR \n",
    "    A = calculate_A(H11)\n",
    "    N2_tilde = ((H12**2)*P1+N2)\n",
    "    NR_tilde = ((H1R**2)*P1+NR)\n",
    "    P_Z = ((H12*H1R*P1)/(np.sqrt(N2_tilde*NR_tilde)))\n",
    "    K1 = ((H2R**2)*N2_tilde+(H22**2)*NR_tilde-2*H2R*H22*P_Z*np.sqrt(N2_tilde*NR_tilde))\n",
    "    K2 = ((1-P_Z**2)*NR_tilde*N2_tilde)\n",
    "\n",
    "\n",
    "    C1 = K1*(HR1**2)*((H22**2)*(HR1**2)-(HR2**2)*(H21**2))\n",
    "    C2 = K1*(HR2**2)*(H21**2)*A-2*(HR1**2)*(H22**2)*K1*A-(HR1**2)*(H22**2)*(H21**2)*K2\n",
    "    C3 = (H22**2)*K1*(A**2)+(H22**2)*(H21**2)*K2*A\n",
    "    C4 = K2*(HR2**2)*(H21**2)**2-N2_tilde*K1*(H21**2)*(HR1**2)\n",
    "    C5 = N2_tilde*(H21**2)*(K1*A+K2*(H21**2))\n",
    "\n",
    "    f = (C1*x**2+C2*x+C3)/(C4*x+C5)\n",
    "    \n",
    "    return f"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Bruteforce"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [],
   "source": [
    "def BF_V2(HR2, H12, H1R, HR1, H22, H2R, H21, H11):\n",
    "    ''' Bruteforce '''\n",
    "    Pr = np.linspace(0.0, PR_max, 1000)\n",
    "    P2 = np.linspace(0.0, P2_max, 1000)\n",
    "\n",
    "    B,C = np.meshgrid(Pr, P2)\n",
    "              \n",
    "    A = calculate_A(H11)     \n",
    "    \n",
    "    # if QoS constraint respected\n",
    "    mask = B <= (A-(H21**2)*C)/(HR1**2)\n",
    "    #mask2 = C <= (A-(HR1**2)*B)/(H21**2)\n",
    "\n",
    "    B = B[mask]\n",
    "    C = C[mask]\n",
    "    \n",
    "    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, B, C)\n",
    "\n",
    "    ind = np.argmax(SNR)\n",
    "\n",
    "    SNR_opt, pr_opt, ps_opt = SNR[ind], B[ind], C[ind]\n",
    "    \n",
    "    #P2_opt = (A-(HR1**2)*pr_opt)/(H21**2)\n",
    "    \n",
    "    #SNR_opt = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, pr_opt, ps_opt)\n",
    "    #SNR_opt = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, pr_opt, ps_opt)\n",
    "\n",
    "\n",
    "    \n",
    "    return pr_opt, ps_opt, SNR_opt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "def f_domain(HR2, H12, H1R, HR1, H22, H2R, H21, H11):\n",
    "    #HR1, H11, H2R, H1R, H22, HR2, H21, H12\n",
    "    A = calculate_A(H11)\n",
    "    if ( (A/(HR1**2)) < PR_max and (A/(H21**2)) < P2_max ) :  # H1\n",
    "        \n",
    "        x_min, x_max = 0.0, A/(HR1**2)        \n",
    "\n",
    "    elif ( (A/(HR1**2)) < PR_max and (A/(H21**2)) > P2_max ):# H2\n",
    "        \n",
    "        x_min, x_max =  (A-(H21**2)*P2_max)/(HR1**2), A/(HR1**2)\n",
    "\n",
    "    elif ( (A/(HR1**2)) > PR_max and (A/(H21**2)) < P2_max ):# H3\n",
    "        \n",
    "        x_min, x_max =  0.0 , PR_max\n",
    "    elif ( (A/(HR1**2)) > PR_max and (A/(H21**2)) > P2_max ) and P2_max*(H21**2) + PR_max*(HR1**2) > A: \n",
    "        \n",
    "        x_min, x_max =  (A-(H21**2)*P2_max)/(HR1**2), PR_max\n",
    "        \n",
    "    return x_min, x_max, A\n",
    "\n",
    "def CF_V4(HR2, H12, H1R, HR1, H22, H2R, H21, H11):\n",
    "    '''\n",
    "    Version Finale\n",
    "    '''\n",
    "        \n",
    "    A = calculate_A(H11)\n",
    "    #x_min, x_max = f_domain(HR2, H12, H1R, HR1, H22, H2R, H21, H11)\n",
    "    if ( (A/(HR1**2)) < PR_max and (A/(H21**2)) < P2_max ) :  # H1\n",
    "\n",
    "        x_min, x_max = 0.0, A/(HR1**2)        \n",
    "\n",
    "    elif ( (A/(HR1**2)) < PR_max and (A/(H21**2)) > P2_max ):# H2\n",
    "\n",
    "        x_min, x_max =  (A-(H21**2)*P2_max)/(HR1**2), A/(HR1**2)\n",
    "\n",
    "    elif ( (A/(HR1**2)) > PR_max and (A/(H21**2)) < P2_max ):# H3\n",
    "\n",
    "        x_min, x_max =  0.0 , PR_max\n",
    "        \n",
    "    elif ( (A/(HR1**2)) > PR_max and (A/(H21**2)) > P2_max ) and P2_max*(H21**2) + PR_max*(HR1**2) > A: #H4\n",
    "\n",
    "        x_min, x_max =  (A-(H21**2)*P2_max)/(HR1**2), PR_max\n",
    "        \n",
    "    C1, C2, C3, C4, C5, delta = Delta(HR2, H12, H1R, HR1, H22, H2R, H21, H11)\n",
    "\n",
    "    try :\n",
    "        if C1*C4>0:\n",
    "            #print('C1*C4>0')\n",
    "            \n",
    "            if delta > 0 :\n",
    "                x_1 = (-(C1*C5)-np.sqrt(delta))/(C1*C4) \n",
    "                x_2 = (-(C1*C5)+np.sqrt(delta))/(C1*C4) \n",
    "                if (x_1 <= x_min <= x_max <= x_2): \n",
    "                    Pr_opt = x_min \n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt,P2_opt)\n",
    "                elif (x_1<=x_min<=x_2<=x_max and f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_min)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_max)): \n",
    "                    Pr_opt = x_min \n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt,P2_opt)\n",
    "                      \n",
    "                elif x_1 >= x_max :\n",
    "                    Pr_opt = x_max \n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "               \n",
    "                elif x_2 <= x_max :\n",
    "                    Pr_opt = x_max \n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "                elif  (x_1 <= x_min <= x_2 <= x_max and f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_max)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_min)):\n",
    "                    Pr_opt = x_max \n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                    \n",
    "                elif (x_min <= x_1 <= x_2<=x_max and f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_max)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_1)): \n",
    "                    Pr_opt = x_max \n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "                elif (x_min <= x_1 <= x_max <=x_2) :\n",
    "                    Pr_opt = x_1\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                elif (x_min<=x_1<=x_2<=x_max and f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_1)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_max)):\n",
    "                    Pr_opt = x_1\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "            elif delta == 0 :\n",
    "                \n",
    "                x_0 = -C5/C4\n",
    "                Pr_opt = x_max\n",
    "                P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                \n",
    "            else :\n",
    "                    Pr_opt = x_max\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                    \n",
    "\n",
    "\n",
    "        elif C1*C4 == 0.0 : \n",
    "            #print('C1*C4 == 0')\n",
    "            x_ = -(C2*C5-C3*C4)/(2*C1*C5)\n",
    "\n",
    "            if (C1*C5)>0.0 :\n",
    "                \n",
    "                if (C2*C5-C3*C4)>0.0 :\n",
    "\n",
    "                    Pr_opt = x_max\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "                else : \n",
    "\n",
    "                    if x_ < x_min : \n",
    "\n",
    "                        Pr_opt = x_max\n",
    "                        P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                        SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "                    elif (x_min <= x_ <= x_max) : \n",
    "\n",
    "                        if f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_min)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_max) :\n",
    "                            \n",
    "                            Pr_opt = x_min\n",
    "                            P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                            SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "                        else : \n",
    "                            \n",
    "                            Pr_opt = x_max\n",
    "                            P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                            SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "                    elif x_ > x_max :\n",
    "\n",
    "                        Pr_opt = x_min\n",
    "                        P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                        SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "            elif (C1*C5) == 0.0 : \n",
    "                #print('C1*C5')\n",
    "                if (C2*C5-C3*C4) > 0.0 :\n",
    "                    \n",
    "                    Pr_opt = x_max\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "                else :\n",
    "                    \n",
    "                    Pr_opt = x_min\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "            else :\n",
    "                if (C2*C5-C3*C4)>0.0 :\n",
    "                    #print((C2*C5-C3*C4))\n",
    "                    if x_ < x_min :\n",
    "                        \n",
    "                        Pr_opt = x_min\n",
    "                        P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                        SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "                    elif (x_min <= x_ <= x_max):\n",
    "                        \n",
    "                        Pr_opt = x_\n",
    "                        P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                        SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "                    elif x_ > x_max :\n",
    "                        \n",
    "                        Pr_opt = x_max  #     x_    \n",
    "                        P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                        SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                \n",
    "                else : # (C2*C5-C3*C4) < 0 ...\n",
    "                        Pr_opt = x_min        \n",
    "                        P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                        SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "        else :\n",
    "            #print('C1*C4<0')\n",
    "\n",
    "            if delta > 0.0 :\n",
    "                #print('delta>0')\n",
    "                x_1 = (-(C1*C5)-np.sqrt(delta))/(C1*C4) \n",
    "                x_2 = (-(C1*C5)+np.sqrt(delta))/(C1*C4)\n",
    "                \n",
    "                if x_1 <= x_min  : \n",
    "                    \n",
    "                    Pr_opt = x_min\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                \n",
    "                elif (x_min<=x_2<=x_max and f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_min)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11,x_max)): \n",
    "                    Pr_opt = x_min\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "                elif  x_max <= x_2:\n",
    "                    Pr_opt = x_min\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                elif  (x_min<=x_2<=x_max<=x_1 and f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_max)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11,x_min)):\n",
    "                    Pr_opt = x_min\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "\n",
    "\n",
    "                elif (x_min <= x_1 <= x_max) :\n",
    "                    \n",
    "                    Pr_opt = x_1\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                elif (x_min <= x_2 <= x_max <= x_1 and f(HR2, H12, H1R, HR1, H22, H2R, H21, H11, x_max)>f(HR2, H12, H1R, HR1, H22, H2R, H21, H11,x_min)) :\n",
    "                    \n",
    "                    Pr_opt = x_max\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                elif (x_2 <= x_min <= x_max <= x_1) :\n",
    "                    \n",
    "                    Pr_opt = x_max\n",
    "                    P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                    SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "               \n",
    "                \n",
    "            elif delta == 0 :\n",
    "                #print('delta=0')\n",
    "\n",
    "\n",
    "                x_0 = -C5/C4\n",
    "                Pr_opt = x_min\n",
    "                P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                \n",
    "            else :\n",
    "                #print('delta<0')\n",
    "                Pr_opt = x_min\n",
    "                P2_opt = (A-(HR1**2)*Pr_opt)/(H21**2)\n",
    "                SNR = f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, Pr_opt, P2_opt)\n",
    "                \n",
    "\n",
    "     \n",
    "\n",
    "    except UnboundLocalError:\n",
    "        \n",
    "        Pr_opt, P2_opt, SNR = PR_max, P2_max ,f_obj(HR2, H12, H1R, HR1, H22, H2R, H21, H11, PR_max, P2_max)\n",
    "        #print('[H5]-[H6]')\n",
    "    \n",
    "    \n",
    "    return Pr_opt, P2_opt, SNR\n",
    "            \n",
    "    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Analytic soluton"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Phase test"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Dataset creation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [],
   "source": [
    "Nbr = 10000\n",
    "pos_min = 0 # min position in cell\n",
    "pos_max = 10 # max position in cell\n",
    "\n",
    "#-----------------------------------------------------------------------------------------------------------#\n",
    "#-----------------------------------------------------------------------------------------------------------#\n",
    "\n",
    "def channel_gain_with_gaussian_fading(d, sigma=7, alpha=3):  # channel gain model\n",
    "  \"\"\"\n",
    "      channel gain model.\n",
    "      Args:\n",
    "         d: distance between source and destination\n",
    "         sigma: gaussian fading sigma\n",
    "         alpha: is the path loss factor\n",
    "      Returns:\n",
    "        channel gain\n",
    "  \"\"\"    \n",
    "  s = np.random.normal(0.0, sigma, d.shape[0])\n",
    "  h = s/np.sqrt(1.0 + np.power(d, alpha))\n",
    "\n",
    "  return h\n",
    "\n",
    "def data_filter(HR2, H12, H1R, HR1, H22, H2R, H21, H11):\n",
    "  s1 = 10**(-40/10)/10\n",
    "  s = np.sqrt(s1)\n",
    "  mask = np.all(np.stack([HR2, H12, H1R, HR1, H22, H2R, H21, H11], axis=1)>=s, axis=1)\n",
    "  return HR2[mask], H12[mask], H1R[mask], HR1[mask], H22[mask], H2R[mask], H21[mask], H11[mask]\n",
    "\n",
    "\n",
    "#-----------------------------------------------------------------------------------------------------------#\n",
    "#-----------------------------------------------------------------------------------------------------------#\n",
    "\n",
    "U_1 = np.random.uniform(pos_min, pos_max, (Nbr, 2))\n",
    "D_1 = np.random.uniform(pos_min, pos_max, (Nbr, 2))\n",
    "U_2 = np.random.uniform(pos_min, pos_max, (Nbr, 2))\n",
    "D_2 = np.random.uniform(pos_min, pos_max, (Nbr, 2))\n",
    "R  = 5 * np.ones((Nbr, 2))\n",
    "\n",
    "d_11 = np.linalg.norm(U_1 - D_1, axis=1)\n",
    "d_12 = np.linalg.norm(U_1 - D_2, axis=1)\n",
    "d_1R = np.linalg.norm(U_1 - R  , axis=1)\n",
    "d_21 = np.linalg.norm(U_2 - D_1, axis=1)\n",
    "d_22 = np.linalg.norm(U_2 - D_2, axis=1)\n",
    "d_2R = np.linalg.norm(U_2 - R  , axis=1)\n",
    "d_R1 = np.linalg.norm(R   - D_1, axis=1)\n",
    "d_R2 = np.linalg.norm(R   - D_2, axis=1)\n",
    "\n",
    "var  = 7\n",
    "\n",
    "h_R2 = channel_gain_with_gaussian_fading(d_R2, var)\n",
    "h_12 = channel_gain_with_gaussian_fading(d_12, var)\n",
    "h_1R = channel_gain_with_gaussian_fading(d_1R, var)\n",
    "h_R1 = channel_gain_with_gaussian_fading(d_R1, var)\n",
    "h_22 = channel_gain_with_gaussian_fading(d_22, var)\n",
    "h_2R = channel_gain_with_gaussian_fading(d_2R, var)\n",
    "h_21 = channel_gain_with_gaussian_fading(d_21, var)\n",
    "h_11 = channel_gain_with_gaussian_fading(d_11, var)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Run CF Code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "CF_res = mult_res_v2(CF_V4, h_R2, h_12, h_1R, h_R1, h_22, h_2R, h_21, h_11)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Run Bruteforce Code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [],
   "source": [
    "BF_res = mult_res_v2(BF_V2, h_R2, h_12, h_1R, h_R1, h_22, h_2R, h_21, h_11)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0"
      ]
     },
     "execution_count": 38,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.size(np.where(BF_res[:, 2] > CF_res[:, 2]))\n",
    "#np.where(BF_res[:, 2] > CF_res[:, 2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Faux !! np.size(np.where(h_21**2*CF_res[:, 1]+h_R1**2*CF_res[:, 0] > calculate_A(h_11)))\n",
    "\n",
    "#CF_V4(h_R2[4339], h_12[4339], h_1R[4339], h_R1[4339], h_22[4339], h_2R[4339], h_21[4339], h_11[4339])\n",
    "#h_21[indx_of_max]**2*CF_res[indx_of_max, 1]+h_R1[indx_of_max]**2*CF_res[indx_of_max, 0] > calculate_A(h_11_n[indx_of_max])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "Pr_distance_max = np.max(np.abs(CF_res[:, 0] - BF_res[:, 0]))\n",
    "P2_distance_max = np.max(np.abs(CF_res[:, 1] - BF_res[:, 1]))\n",
    "SNR_distance_max = np.max(np.abs(CF_res[:, 2] - BF_res[:, 2]))\n",
    "distance_SNR = np.abs(CF_res[:, 2] - BF_res[:, 2])\n",
    "#distance_SNR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.0154347411735927"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "indx_of_max = np.where(np.max(distance_SNR)==distance_SNR)\n",
    "np.max(distance_SNR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PR_SA : [[0.]]\n",
      "PR_BF : [[0.06006006]]\n",
      "P2_SA : [[0.26685174]]\n",
      "P2_BF : [[0.26026026]]\n",
      "SNR_SA : [[41.12131901]]\n",
      "SNR_BF : [[40.10588427]]\n"
     ]
    }
   ],
   "source": [
    "print(\"PR_SA :\",CF_res[indx_of_max, 0])\n",
    "print(\"PR_BF :\",BF_res[indx_of_max, 0])\n",
    "print(\"P2_SA :\",CF_res[indx_of_max, 1])\n",
    "print(\"P2_BF :\",BF_res[indx_of_max, 1])\n",
    "print(\"SNR_SA :\",CF_res[indx_of_max, 2])\n",
    "print(\"SNR_BF :\",BF_res[indx_of_max, 2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0.12711753]])"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# changement de l'echelle de la bruteforce pour verifier ce cas (1000 a 10000)\n",
    "S = BF_V2(h_R2[indx_of_max], h_12[indx_of_max], h_1R[indx_of_max], h_R1[indx_of_max], h_22[indx_of_max], h_2R[indx_of_max], h_21[indx_of_max], h_11[indx_of_max])\n",
    "# recalcul de la distance\n",
    "np.abs(S[2] - CF_res[indx_of_max, 2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Analytic Solution :\n",
      "   PR    P2   SNR\n",
      "[[10.  , 5.15, 1.02],\n",
      " [ 0.  , 1.77, 0.22],\n",
      " [ 0.  , 5.03, 0.16],\n",
      " ...,\n",
      " [10.  , 5.07,14.01],\n",
      " [10.  , 0.08, 0.  ],\n",
      " [10.  , 6.3 , 1.16]]\n",
      "--------------------------------------------------\n",
      "Brute Force : \n",
      "   PR    P2   SNR\n",
      "[[10.  , 5.15, 1.02],\n",
      " [ 0.06, 1.76, 0.22],\n",
      " [ 0.  , 5.03, 0.16],\n",
      " ...,\n",
      " [10.  , 5.07,13.99],\n",
      " [10.  , 0.07, 0.  ],\n",
      " [10.  , 6.29, 1.16]]\n"
     ]
    }
   ],
   "source": [
    "#np.set_pxrintoptions(formatter={'float': lambda x : '{0:0.2f}'.format(x)})\n",
    "print('Analytic Solution :')\n",
    "print('   PR    P2   SNR')\n",
    "print(np.array2string(CF_res, precision=2, separator=',',\n",
    "                       suppress_small=True))\n",
    "\n",
    "print('--------------------------------------------------')\n",
    "\n",
    "print('Brute Force : ')\n",
    "print('   PR    P2   SNR')\n",
    "print(np.array2string(BF_res, precision=2, separator=',', suppress_small=True))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 47,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "10.000000000485224"
      ]
     },
     "execution_count": 47,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "Pr_distance_max = np.max(np.abs(CF_res[:, 0] - BF_res[:, 0]))\n",
    "P2_distance_max = np.max(np.abs(CF_res[:, 1] - BF_res[:, 1]))\n",
    "SNR_distance_max = np.max(np.abs(CF_res[:, 2] - BF_res[:, 2]))\n",
    "\n",
    "np.max(CF_res[:, 1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Pr</th>\n",
       "      <th>P2</th>\n",
       "      <th>SNR</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>10000.000000</td>\n",
       "      <td>10000.000000</td>\n",
       "      <td>1.000000e+04</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>5.279746</td>\n",
       "      <td>7.471945</td>\n",
       "      <td>8.597277e+00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>4.318922</td>\n",
       "      <td>3.714079</td>\n",
       "      <td>5.538911e+01</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.006965</td>\n",
       "      <td>1.525856e-07</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>0.569573</td>\n",
       "      <td>4.124612</td>\n",
       "      <td>8.799407e-02</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>4.878257</td>\n",
       "      <td>10.000000</td>\n",
       "      <td>5.491686e-01</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>10.000000</td>\n",
       "      <td>10.000000</td>\n",
       "      <td>2.572905e+00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>10.000000</td>\n",
       "      <td>10.000000</td>\n",
       "      <td>1.796350e+03</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 Pr            P2           SNR\n",
       "count  10000.000000  10000.000000  1.000000e+04\n",
       "mean       5.279746      7.471945  8.597277e+00\n",
       "std        4.318922      3.714079  5.538911e+01\n",
       "min        0.000000      0.006965  1.525856e-07\n",
       "25%        0.569573      4.124612  8.799407e-02\n",
       "50%        4.878257     10.000000  5.491686e-01\n",
       "75%       10.000000     10.000000  2.572905e+00\n",
       "max       10.000000     10.000000  1.796350e+03"
      ]
     },
     "execution_count": 45,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "CF_res_pd = pd.DataFrame({'Pr': CF_res[:, 0], 'P2': CF_res[:, 1], 'SNR': CF_res[:, 2] })\n",
    "CF_res_pd.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Pr</th>\n",
       "      <th>P2</th>\n",
       "      <th>SNR</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>count</th>\n",
       "      <td>10000.000000</td>\n",
       "      <td>10000.000000</td>\n",
       "      <td>1.000000e+04</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>mean</th>\n",
       "      <td>5.297651</td>\n",
       "      <td>7.469982</td>\n",
       "      <td>8.595267e+00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>std</th>\n",
       "      <td>4.302517</td>\n",
       "      <td>3.714411</td>\n",
       "      <td>5.538563e+01</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>min</th>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.010010</td>\n",
       "      <td>1.521265e-07</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25%</th>\n",
       "      <td>0.620621</td>\n",
       "      <td>4.134134</td>\n",
       "      <td>8.790783e-02</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50%</th>\n",
       "      <td>4.934935</td>\n",
       "      <td>10.000000</td>\n",
       "      <td>5.488385e-01</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>75%</th>\n",
       "      <td>10.000000</td>\n",
       "      <td>10.000000</td>\n",
       "      <td>2.571248e+00</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>max</th>\n",
       "      <td>10.000000</td>\n",
       "      <td>10.000000</td>\n",
       "      <td>1.796350e+03</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                 Pr            P2           SNR\n",
       "count  10000.000000  10000.000000  1.000000e+04\n",
       "mean       5.297651      7.469982  8.595267e+00\n",
       "std        4.302517      3.714411  5.538563e+01\n",
       "min        0.000000      0.010010  1.521265e-07\n",
       "25%        0.620621      4.134134  8.790783e-02\n",
       "50%        4.934935     10.000000  5.488385e-01\n",
       "75%       10.000000     10.000000  2.571248e+00\n",
       "max       10.000000     10.000000  1.796350e+03"
      ]
     },
     "execution_count": 46,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "\n",
    "BF_res_pd = pd.DataFrame({'Pr': BF_res[:, 0], 'P2': BF_res[:, 1], 'SNR': BF_res[:, 2] })\n",
    "BF_res_pd.describe()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 132,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import mean_absolute_error, mean_squared_error\n",
    "\n",
    "y_true, y_predicted = CF_res[:,2].astype(\"float\"), BF_res[:,2].astype(\"float\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  Mean Absolute Error (SNR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 133,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.003648147975755582"
      ]
     },
     "execution_count": 133,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mean_absolute_error(y_true, y_predicted)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##  Mean Squared Error (SNR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 134,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.00471967919039298"
      ]
     },
     "execution_count": 134,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mean_squared_error(y_true, y_predicted)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.8.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}