P_cg_gr_bt_conv_v32.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "import  pandas               as     pd\n",
    "import  pandapower           as     pp\n",
    "import  numpy                as     np\n",
    "import  math                 as     mt\n",
    "import  cmath \n",
    "import  matplotlib.pyplot    as     plt\n",
    "import  time                 as     tm\n",
    "from    scipy.interpolate    import griddata\n",
    "from    scipy.optimize       import minimize\n",
    "from    scipy.optimize       import least_squares\n",
    "from    random               import uniform\n",
    "from    copy                 import deepcopy\n",
    "from    itertools            import chain"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set a Element Class that is common to all the elements in the grid\n",
    "\n",
    "class Node():\n",
    "    def __init__(self, node_number, node_name, vn_kv, elec_type):\n",
    "        self.node_number = node_number\n",
    "        self.node_name = node_name\n",
    "        self.vn_kv = vn_kv\n",
    "        self.elec_type = elec_type\n",
    "\n",
    "class Line():\n",
    "    def __init__(self, from_bus, to_bus, line_id, l_km, r_ohm_km, x_ohm_km, c_ohm_km, Imax_kA):\n",
    "        self.from_bus = from_bus\n",
    "        self.to_bus = to_bus\n",
    "        self.line_id = line_id\n",
    "        self.l_km = l_km\n",
    "        self.r_ohm_km = r_ohm_km\n",
    "        self.x_ohm_km = x_ohm_km\n",
    "        self.c_ohm_km = c_ohm_km\n",
    "        self.Imax_kA = Imax_kA\n",
    "        \n",
    "class Transformer():\n",
    "    def __init__(self, from_bus, to_bus, trafo_id, sn, vk_percent, vkr_percent, pfe_kw, i0_percent, shift_degree, vector_group, tap_side, min_tap, max_tap, tap_step_degree, tap_step, tap_ps, vk0_percent, vkr0_percent, mag0_percent, mag0_rx, si0_hv_partial):\n",
    "        self.from_bus = from_bus\n",
    "        self.to_bus = to_bus\n",
    "        self.trafo_id = trafo_id\n",
    "        self.sn = sn\n",
    "        self.vk_percent = vk_percent\n",
    "        self.vkr_percent = vkr_percent\n",
    "        self.pfe_kw = pfe_kw\n",
    "        self.i0_percent = i0_percent\n",
    "        self.shift_degree = shift_degree\n",
    "        self.vector_group = vector_group\n",
    "        self.tap_side = tap_side\n",
    "        self.min_tap = min_tap\n",
    "        self.max_tap = max_tap\n",
    "        self.tap_step_degree = tap_step_degree\n",
    "        self.tap_step = tap_step\n",
    "        self.tap_ps = tap_ps\n",
    "        self.vk0_percent = vk_percent\n",
    "        self.vkr0_percent = vkr_percent\n",
    "        self.mag0_percent = i0_percent\n",
    "        self.mag0_rx = mag0_rx\n",
    "        self.si0_hv_partial = si0_hv_partial\n",
    "\n",
    "class Converter():\n",
    "    def __init__(self, from_bus, to_bus, conv_id, convtype, pn, r_ohm, x_ohm, r_dc_ohm, pl_dc_mw, control_mode, control_value, eff_table):\n",
    "        self.from_bus = from_bus\n",
    "        self.to_bus = to_bus\n",
    "        self.conv_id = conv_id\n",
    "        self.convtype = convtype\n",
    "        self.pn = pn\n",
    "        self.r_ohm = r_ohm\n",
    "        self.x_ohm = x_ohm\n",
    "        self.r_dc_ohm = r_dc_ohm\n",
    "        self.pl_dc_mw = pl_dc_mw\n",
    "        self.control_mode = control_mode\n",
    "        self.control_value = control_value\n",
    "        self.eff_table = eff_table\n",
    "        \n",
    "\n",
    "class Elem:\n",
    "    def __init__(self, owner, bus):\n",
    "        self.owner = owner\n",
    "        self.bus = bus\n",
    "\n",
    "class Load(Elem):\n",
    "    def __init__(self, owner, bus, load_id, ctype, dp_cont, vn, pn, droop):\n",
    "        super().__init__(owner, bus)\n",
    "        self.load_id = load_id\n",
    "        self.ctype = ctype\n",
    "        self.dp_cont = dp_cont\n",
    "        self.vn = vn\n",
    "        self.pn = pn\n",
    "        self.droop = droop\n",
    "\n",
    "\n",
    "class Gen(Elem):\n",
    "    def __init__(self, owner, bus, gen_id, gentype, dp_cont, vn, pn, droop):\n",
    "        super().__init__(owner, bus)\n",
    "        self.gen_id = gen_id\n",
    "        self.gentype = gentype\n",
    "        self.dp_cont = dp_cont\n",
    "        self.vn = vn\n",
    "        self.pn = pn\n",
    "        self.droop = droop\n",
    "\n",
    "class Stor(Elem):\n",
    "    def __init__(self, owner, bus, stor_id, stortype, dp_cont, ini_SOC, vn, pn, max_ener, max_p_discharge, max_p_charge, max_SOC, min_SOC, droop):\n",
    "        super().__init__(owner, bus)\n",
    "        self.stor_id = stor_id\n",
    "        self.stortype = stortype\n",
    "        self.dp_cont = dp_cont\n",
    "        self.ini_SOC = ini_SOC\n",
    "        self.vn = vn\n",
    "        self.pn = pn\n",
    "        self.max_ener = max_ener\n",
    "        self.max_p_discharge = max_p_discharge\n",
    "        self.max_p_charge = max_p_charge\n",
    "        self.max_SOC = max_SOC\n",
    "        self.min_SOC = min_SOC\n",
    "        self.droop = droop\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Read the Excel file\n",
    "print(\">Read the Excel file\")\n",
    "xl_file = pd.ExcelFile('C:/Users/HP/Desktop/Engenharia Eletrotécnica/Doutoramento/1º Ano/Tese/Modelação/Redes de teste/4_grid_v6.xlsx')\n",
    "# xl_file = pd.ExcelFile('grid_data_buildingtest.xlsx')\n",
    "\n",
    "# Load data from Excel sheets\n",
    "print(\">Load data from Excel sheets\")\n",
    "general_data = xl_file.parse('GeneralInfo')\n",
    "node_data = xl_file.parse('NodeInfo')\n",
    "line_data = xl_file.parse('LineInfo')\n",
    "transformer_data = xl_file.parse('TransformerInfo')\n",
    "converter_data = xl_file.parse('ConverterInfo')\n",
    "load_data = xl_file.parse('LoadInfo',header = 0)#, header = [0,1,2],index_col=[0,2]) \n",
    "gen_data = xl_file.parse('GenInfo', header = 0) \n",
    "stor_data = xl_file.parse('StorInfo', header = 0) \n",
    "\n",
    "#print(general_data)\n",
    "#print(node_data)\n",
    "#print(line_data)\n",
    "#print(transformer_data)\n",
    "#print(converter_data)\n",
    "#print(load_data)\n",
    "#print(gen_data)\n",
    "#print(stor_data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Grid data for time variant PF -> Creation of the main tables of network data\n",
    "\n",
    "def Griddata(general_data,node_data,line_data,transformer_data,converter_data,load_data,gen_data,stor_data):\n",
    "\n",
    "    # Ensure the column names match exactly and there is no leading/trailing whitespace\n",
    "    general_data.columns = general_data.columns.str.strip()\n",
    "    # Define the base voltage and power from the baseValue sheet\n",
    "    #print(\">Define the base voltage and power from the baseValue sheet\")\n",
    "    sim_periods = general_data[general_data['General Information'].str.strip() == 'Simulation Periods']['Value']\n",
    "    period_duration = general_data[general_data['General Information'].str.strip() == 'Periods Duration (h)']['Value']\n",
    "    base_power = general_data[general_data['General Information'].str.strip() == 'Base Power (kW)']['Value']\n",
    "    slack_comp = general_data[general_data['General Information'].str.strip() == 'Slack Gen Component']['Value']\n",
    "\n",
    "    if sim_periods.empty:\n",
    "        raise ValueError(\"The value of Number of Periods is missing.\")\n",
    "    elif period_duration.empty:\n",
    "        raise ValueError(\"The value of Period Duration is missing.\")\n",
    "    elif base_power.empty:\n",
    "        raise ValueError(\"Base value for Power is missing.\")\n",
    "    elif slack_comp.empty:\n",
    "        raise ValueError(\"The value of Slack Component is missing.\")\n",
    "\n",
    "    base_power_kW = base_power.values[0]\n",
    "    #print(base_power)\n",
    "\n",
    "    ## Buses\n",
    "\n",
    "    Nodelist = [Node(0,i,0,0) for i in node_data['Name']]\n",
    "\n",
    "    for i in range(0,len(Nodelist)): \n",
    "        \n",
    "        j = node_data.index[Nodelist[i].node_name == node_data['Name']]\n",
    "\n",
    "        Nodelist[i].node_number = node_data['Node Number'].iloc[j]\n",
    "    \n",
    "        Nodelist[i].node_name = node_data['Name'].iloc[j]       \n",
    "    \n",
    "        Nodelist[i].vn_kv = node_data['Vn (kV)'].iloc[j]\n",
    "    \n",
    "        Nodelist[i].elec_type = node_data['Type'].iloc[j]\n",
    "\n",
    "    ## Network equipment\n",
    "    #print(\">Display network data\")\n",
    "    #print(network_data)\n",
    "\n",
    "    # Line equipment\n",
    "\n",
    "    Linelist = [Line(0, 0, i, 0, 0, 0, 0, 0) for i in line_data['Component']]\n",
    "\n",
    "    for i in range(0,len(Linelist)): \n",
    "        \n",
    "        j = line_data.index[Linelist[i].line_id == line_data['Component']]\n",
    "\n",
    "        Linelist[i].from_bus = line_data['Node_i'].iloc[j]\n",
    "    \n",
    "        Linelist[i].to_bus = line_data['Node_j'].iloc[j]       \n",
    "    \n",
    "        Linelist[i].l_km = line_data['Length_ij_km'].iloc[j]\n",
    "    \n",
    "        Linelist[i].r_ohm_km = line_data['Rij_ohm_km'].iloc[j]\n",
    "    \n",
    "        Linelist[i].x_ohm_km = line_data['Xij_ohm_km'].iloc[j]                     \n",
    "    \n",
    "        Linelist[i].c_ohm_km = line_data['Cij_ohm_km'].iloc[j]\n",
    "    \n",
    "        Linelist[i].Imax_kA = line_data['Imax_ij_kA'].iloc[j] \n",
    "        \n",
    "    # Transformer equipment\n",
    "    \n",
    "    Trafolist = [Transformer(0, 0, i, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0) for i in transformer_data['Component']]\n",
    "\n",
    "    for i in range(0,len(Trafolist)): # This loops assume that the Load ID is coherent with the index its position (Miss make the data given to each load being dependent of its ID instead of its position in table)\n",
    "        \n",
    "        j = transformer_data.index[Trafolist[i].trafo_id == transformer_data['Component']]\n",
    "\n",
    "        Trafolist[i].from_bus = transformer_data['Node_i'].iloc[j]\n",
    "    \n",
    "        Trafolist[i].to_bus = transformer_data['Node_j'].iloc[j]       \n",
    "    \n",
    "        Trafolist[i].sn = transformer_data['Sn (MVA)'].iloc[j]\n",
    "    \n",
    "        Trafolist[i].vk_percent = transformer_data['Short-circuit voltage (pu)'].iloc[j]\n",
    "    \n",
    "        Trafolist[i].vkr_percent = transformer_data['Resistive part of Short-circuit voltage (pu)'].iloc[j]\n",
    "        \n",
    "        Trafolist[i].pfe_kw = transformer_data['Iron losses (kW)'].iloc[j]\n",
    "    \n",
    "        Trafolist[i].i0_percent = transformer_data['Magnetization Current (pu)'].iloc[j]\n",
    "    \n",
    "        Trafolist[i].shift_degree = transformer_data['Shift degree'].iloc[j]\n",
    "\n",
    "        Trafolist[i].vector_group = transformer_data['Vector Group'].iloc[j]\n",
    "\n",
    "        Trafolist[i].tap_side = transformer_data['Tap Side'].iloc[j]      \n",
    "\n",
    "        Trafolist[i].min_tap = transformer_data['Minimal Tap'].iloc[j]      \n",
    "                \n",
    "        Trafolist[i].max_tap = transformer_data['Maximum Tap'].iloc[j]      \n",
    "\n",
    "        Trafolist[i].tap_step_degree = transformer_data['Tap step degree'].iloc[j]\n",
    "\n",
    "        Trafolist[i].tap_step = transformer_data['Tap step (%)'].iloc[j]\n",
    "\n",
    "        Trafolist[i].tap_ps = transformer_data['Tap Phase Shifter'].iloc[j]\n",
    "    \n",
    "        Trafolist[i].vk0_percent = transformer_data['Zero sequence short-circuit voltage (pu)'].iloc[j]    \n",
    "\n",
    "        Trafolist[i].vk0r_percent = transformer_data['Zero sequence resistive part of short-circuit voltage (pu)'].iloc[j]     \n",
    "        \n",
    "        Trafolist[i].mag0_percent = transformer_data['Zero sequence of magnetization current (pu)'].iloc[j]  \n",
    "\n",
    "        Trafolist[i].mag0_rx = transformer_data['Ratio R/X of Zero sequence of Magnetization Current'].iloc[j]    \n",
    "\n",
    "        Trafolist[i].si0_hv_partial = transformer_data['Zero sequence short-circuit impedance distribution in HV side'].iloc[j]     \n",
    "           \n",
    "    # Converter equipment (In progress)\n",
    "\n",
    "    convCarac = converter_data[['Converter ID','Unnamed: 1','Unnamed: 2','Unnamed: 3']]\n",
    "    convCarac = convCarac.rename(columns={'Converter ID':convCarac['Converter ID'][0],'Unnamed: 1':convCarac['Unnamed: 1'][0],'Unnamed: 2':convCarac['Unnamed: 2'][0],'Unnamed: 3':convCarac['Unnamed: 3'][0]})\n",
    "    convCarac = convCarac[1:]\n",
    "    convCarac.index = range(0,len(convCarac['Component']))\n",
    "    aux_idx = []\n",
    "    nConvs = max(convCarac['Component'].unique())\n",
    "    for i,_ in enumerate(convCarac.iterrows()):\n",
    "        if (pd.isna(convCarac.iloc[i]).all()):\n",
    "            aux_idx.append(i)\n",
    "    convCarac = convCarac.drop(aux_idx)\n",
    "    convCarac.index = range(0,convCarac.shape[0])\n",
    "    nConvCarac = pd.DataFrame([convCarac['Characteristics ID'].iloc[idx] for idx,_ in enumerate(convCarac['Characteristics ID']) if idx+1 not in convCarac['Characteristics ID'] or type(convCarac['Component'].iloc[idx+1]) is int],index=convCarac['Component'].dropna().unique(),columns=['Number Caracteristics'])\n",
    "    Convlist = [Converter(0,0,i,0,0,0,0,0,0,0,0,0) for i in nConvCarac.index]\n",
    "    aux_nConvCarac = pd.DataFrame(np.cumsum(nConvCarac.values),index=nConvCarac.index)\n",
    "    pos_Comp = np.asarray(np.where(converter_data == 'Component'))\n",
    "    pos_EffCurve = np.asarray(np.where(converter_data == 'Efficiency Curve'))\n",
    "    \n",
    "    for i in nConvCarac.index: # This loops assume that the Converter ID is coherent with the index its position (Miss make the data given to each load being dependent of its ID instead of its position in table)\n",
    "        if i-1 not in nConvCarac.index:\n",
    "            j = 0\n",
    "            jlim = (aux_nConvCarac.loc[i].values-1).item()\n",
    "            while j<=jlim:\n",
    "                if convCarac['Characteristics D.'].iloc[j]=='From Bus':\n",
    "                    Convlist[i-1].from_bus = convCarac['Characteristics V.'].iloc[j]\n",
    "                elif convCarac['Characteristics D.'].iloc[j]=='To Bus':\n",
    "                    Convlist[i-1].to_bus = convCarac['Characteristics V.'].iloc[j]\n",
    "                elif convCarac['Characteristics D.'].iloc[j]=='Converter Type':                  \n",
    "                    Convlist[i-1].convtype = convCarac['Characteristics V.'].iloc[j]    \n",
    "                    j_aux = 0\n",
    "                    if Convlist[i-1].convtype == 'Converter AC/DC':\n",
    "                        while j_aux<=jlim:\n",
    "                            if convCarac['Characteristics D.'].iloc[j_aux]=='Transformer Resistance (Ohm)':\n",
    "                                Convlist[i-1].r_ohm = convCarac['Characteristics V.'].iloc[j_aux]\n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='Transformer Reactance (Ohm)':\n",
    "                                Convlist[i-1].x_ohm = convCarac['Characteristics V.'].iloc[j_aux]\n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='No-Load DC Power Losses (kW)':                  \n",
    "                                Convlist[i-1].pl_dc_mw = convCarac['Characteristics V.'].iloc[j_aux]    \n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='Control Mode':                  \n",
    "                                Convlist[i-1].control_mode = convCarac['Characteristics V.'].iloc[j_aux] \n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='Control Value':                  \n",
    "                                Convlist[i-1].control_value = convCarac['Characteristics V.'].iloc[j_aux] \n",
    "                            j_aux+=1\n",
    "                elif convCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Convlist[i-1].pn = convCarac['Characteristics V.'].iloc[j]\n",
    "                j+=1\n",
    "        else:\n",
    "            j = (aux_nConvCarac.loc[i-1].values).item()\n",
    "            jlim = (aux_nConvCarac.loc[i].values-1).item()\n",
    "            while j<=jlim:\n",
    "                if convCarac['Characteristics D.'].iloc[j]=='From Bus':\n",
    "                    Convlist[i-1].from_bus = convCarac['Characteristics V.'].iloc[j]\n",
    "                elif convCarac['Characteristics D.'].iloc[j]=='To Bus':\n",
    "                    Convlist[i-1].to_bus = convCarac['Characteristics V.'].iloc[j]\n",
    "                elif convCarac['Characteristics D.'].iloc[j]=='Converter Type':                  \n",
    "                    Convlist[i-1].convtype = convCarac['Characteristics V.'].iloc[j]    \n",
    "                    j_aux = 0\n",
    "                    if Convlist[i-1].convtype == 'Converter AC/DC':\n",
    "                        while j_aux<=jlim:\n",
    "                            if convCarac['Characteristics D.'].iloc[j_aux]=='Transformer Resistance (Ohm)':\n",
    "                                Convlist[i-1].r_ohm = convCarac['Characteristics V.'].iloc[j_aux]\n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='Transformer Reactance (Ohm)':\n",
    "                                Convlist[i-1].x_ohm = convCarac['Characteristics V.'].iloc[j_aux]\n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='No-Load DC Power Losses (kW)':                  \n",
    "                                Convlist[i-1].pl_dc_mw = convCarac['Characteristics V.'].iloc[j_aux]    \n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='Control Mode':                  \n",
    "                                Convlist[i-1].control_mode = convCarac['Characteristics V.'].iloc[j_aux] \n",
    "                            elif convCarac['Characteristics D.'].iloc[j_aux]=='Control Value':                  \n",
    "                                Convlist[i-1].control_value = convCarac['Characteristics V.'].iloc[j_aux] \n",
    "                            j_aux+=1\n",
    "                elif convCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Convlist[i-1].pn = convCarac['Characteristics V.'].iloc[j]\n",
    "                j+=1\n",
    "        \n",
    "    # Get the efficency map of each converter to the list Convlist to perform the interpolations in control steps\n",
    "    convEff = converter_data[converter_data.columns[int(pos_EffCurve[1]):]]\n",
    "    convEff.insert(0,converter_data.columns[int(pos_Comp[1])],converter_data[converter_data.columns[int(pos_Comp[1])]])\n",
    "    convEff.columns = convEff.iloc[0,:]   \n",
    "    convEff = convEff[1:]  \n",
    "    convEff.index = range(0,convEff.shape[0])\n",
    "    aux_PF = convEff['Efficiency Curve'].dropna().unique()\n",
    "\n",
    "    for i in range(0,len(Convlist)):\n",
    "        ji = convEff.index[convEff['Component']==(i+1)].item()\n",
    "        jf = ji + len(aux_PF)\n",
    "        Convlist[i].eff_table = convEff.iloc[ji:jf,1:].dropna()        \n",
    "    \n",
    "    print(convCarac)\n",
    "\n",
    "    # Load equipment\n",
    "    #print(\">Display load data\")\n",
    "    #print(load_data)\n",
    "\n",
    "    # Get the caracteristics of the load components in grid from load_data DataFrame\n",
    "\n",
    "    loadCarac = load_data[['Load ID','Unnamed: 1','Unnamed: 2','Unnamed: 3']]\n",
    "    loadCarac = loadCarac.rename(columns={'Load ID':loadCarac['Load ID'][1],'Unnamed: 1':loadCarac['Unnamed: 1'][1],'Unnamed: 2':loadCarac['Unnamed: 2'][1],'Unnamed: 3':loadCarac['Unnamed: 3'][1]})\n",
    "    loadCarac = loadCarac[2:]\n",
    "    loadCarac.index = range(0,len(loadCarac['Component']))\n",
    "    aux_idx = []\n",
    "    nLoads = max(loadCarac['Component'].unique())\n",
    "    for i,_ in enumerate(loadCarac.iterrows()):\n",
    "        if (pd.isna(loadCarac.iloc[i]).all()):\n",
    "            aux_idx.append(i)\n",
    "    loadCarac = loadCarac.drop(aux_idx)\n",
    "    loadCarac.index = range(0,loadCarac.shape[0])\n",
    "    nLoadCarac = pd.DataFrame([loadCarac['Characteristics ID'].iloc[idx] for idx,_ in enumerate(loadCarac['Characteristics ID']) if idx+1 not in loadCarac['Characteristics ID'] or type(loadCarac['Component'].iloc[idx+1]) is int],index=loadCarac['Component'].dropna().unique(),columns=['Number Caracteristics'])\n",
    "    Loadlist = [Load(0,0,i,0,0,0,0,0) for i in nLoadCarac.index]\n",
    "    aux_nLoadCarac = pd.DataFrame(np.cumsum(nLoadCarac.values),index=nLoadCarac.index)\n",
    "\n",
    "    for i in nLoadCarac.index: # This loops assume that the Load ID is coherent with the index its position (Miss make the data given to each load being dependent of its ID instead of its position in table)\n",
    "        if i-1 not in nLoadCarac.index:\n",
    "            j = 0\n",
    "            jlim = (aux_nLoadCarac.loc[i].values-1).item()\n",
    "            while j<=jlim:\n",
    "                if loadCarac['Characteristics D.'].iloc[j]=='Bus Location':\n",
    "                    Loadlist[i-1].bus = loadCarac['Characteristics V.'].iloc[j]\n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Droop Contribution':\n",
    "                    Loadlist[i-1].dp_cont = loadCarac['Characteristics V.'].iloc[j]\n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Charge Type':                  \n",
    "                    Loadlist[i-1].ctype = loadCarac['Characteristics V.'].iloc[j]                      \n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Nominal Voltage (V)':\n",
    "                    Loadlist[i-1].vn = loadCarac['Characteristics V.'].iloc[j]\n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Loadlist[i-1].pn = loadCarac['Characteristics V.'].iloc[j]\n",
    "                #    if loadCarac['Characteristics V.'].iloc[j]=='Car Charger': # Misses the verification of the specifications of specific type of loads \n",
    "                #        auxMinCurrent = 8 # Verify how to read automatically the value\n",
    "                j+=1\n",
    "        else:\n",
    "            j = (aux_nLoadCarac.loc[i-1].values).item()\n",
    "            jlim = (aux_nLoadCarac.loc[i].values-1).item()\n",
    "            while j<=jlim:\n",
    "                if loadCarac['Characteristics D.'].iloc[j]=='Bus Location':\n",
    "                    Loadlist[i-1].bus = loadCarac['Characteristics V.'].iloc[j]\n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Droop Contribution':\n",
    "                    Loadlist[i-1].dp_cont = loadCarac['Characteristics V.'].iloc[j]\n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Charge Type':                   \n",
    "                    Loadlist[i-1].ctype = loadCarac['Characteristics V.'].iloc[j]                     \n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Nominal Voltage (V)':\n",
    "                    Loadlist[i-1].vn = loadCarac['Characteristics V.'].iloc[j]\n",
    "                elif loadCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Loadlist[i-1].pn = loadCarac['Characteristics V.'].iloc[j]\n",
    "                j+=1\n",
    "\n",
    "\n",
    "    #print(loadCarac)\n",
    "\n",
    "    # Get a matrix of all voltages and powers points of Current/OS voltage/power droop curves from load_data DataFrame to perform the interpolations\n",
    "\n",
    "    aux_namCol = ['Load ID','Unnamed: 4','Total Time (h)']\n",
    "    for i,j in enumerate(load_data.columns):\n",
    "        if type(j) is int:                       \n",
    "            aux_namCol.append(j)                    # Gives the correct column names to search the droop information in time in load_data DataFrame \n",
    "\n",
    "\n",
    "    loadDroop = load_data[aux_namCol] \n",
    "    loadDroop = loadDroop.rename(columns={'Load ID':loadDroop['Load ID'][1],'Unnamed: 4':loadDroop['Unnamed: 4'][1],'Total Time (h)':loadDroop['Total Time (h)'][1]})\n",
    "    loadDroop = loadDroop[2:]\n",
    "    maxDroopsLoad = max(loadDroop['Droop ID'].unique())\n",
    "    loadDroop.index = range(0,loadDroop.shape[0])\n",
    "\n",
    "    #print(loadDroop)\n",
    "\n",
    "\n",
    "    for i in nLoadCarac.index:\n",
    "        aux_indx = (np.array(loadDroop.index[loadDroop['Component']==i].tolist())).item()\n",
    "        if aux_indx is not None:\n",
    "            for j in range(0,len(Loadlist)):\n",
    "                if Loadlist[j].load_id == i:\n",
    "                    Loadlist[j].droop = loadDroop.loc[aux_indx:aux_indx+maxDroopsLoad-1]\n",
    "\n",
    "    # This automatic verification of the name pairs of V,P allows to use several numerations of the voltage droop control curves, for exemple, (V1, P1) or (Va, Pa) or (VA,PA)\n",
    "\n",
    "    for i in range(0,len(Loadlist)):                # Sorting process (load by load)\n",
    "        \n",
    "        auxV = pd.DataFrame(Loadlist[i].droop.loc[Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains('V')],'Droop'],index=Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains('V')].tolist())      # List of the names of the voltage points of the Current/OS voltage droop control curves of each load\n",
    "        auxV.sort_values(by='Droop')                                                                                                                                                                                        # This command guarrantees that the voltage names are sorted\n",
    "\n",
    "        auxP = pd.DataFrame(Loadlist[i].droop.loc[Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)],'Droop'],index=Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)].tolist())      # List of the names of the voltage points of the Current/OS voltage droop control curves of each load\n",
    "        auxP.sort_values(by='Droop')                                                                                                                                                                                        # This command guarrantees that the voltage names are sorted\n",
    "\n",
    "        auxVP = list(zip(auxV.index,auxP.index))    # List with tuples with the indexes corresponding to the names of the \"coordinates\" of the points of the Current/OS voltage droop control curves of each load\n",
    "\n",
    "        for k in Loadlist[i].droop.columns:\n",
    "            if type(k) is int:\n",
    "                droop_list = []\n",
    "                for j,_ in enumerate(auxVP): \n",
    "                    droop_list.append(Loadlist[i].droop.loc[np.array(list(auxVP[j])),k].values)  \n",
    "                df = pd.DataFrame(droop_list, columns=['x','y'])\n",
    "                df_sorted = df.sort_values(by=['x','y'])\n",
    "                Loadlist[i].droop.loc[auxV.index,k] = df_sorted['x'].values\n",
    "                Loadlist[i].droop.loc[auxP.index,k] = df_sorted['y'].values\n",
    "\n",
    "    # Set the Load Class to characterize each load\n",
    "\n",
    "    # Generator equipment\n",
    "    #print(\">Display generator data\")\n",
    "    #print(gen_data)\n",
    "\n",
    "    # Get the caracteristics of the generator components in grid from gen_data DataFrame\n",
    "\n",
    "    genCarac = gen_data[['Gen ID','Unnamed: 1','Unnamed: 2','Unnamed: 3']]\n",
    "    genCarac = genCarac.rename(columns={'Gen ID':genCarac['Gen ID'][1],'Unnamed: 1':genCarac['Unnamed: 1'][1],'Unnamed: 2':genCarac['Unnamed: 2'][1],'Unnamed: 3':genCarac['Unnamed: 3'][1]})\n",
    "    genCarac = genCarac[2:]\n",
    "    genCarac.index = range(0,len(genCarac['Component']))\n",
    "    aux_idx = []\n",
    "    nGens = max(genCarac['Component'].unique())\n",
    "    for i,_ in enumerate(genCarac.iterrows()):\n",
    "        if (pd.isna(genCarac.iloc[i]).all()):\n",
    "            aux_idx.append(i)\n",
    "    genCarac = genCarac.drop(aux_idx)\n",
    "    genCarac.index = range(0,genCarac.shape[0])\n",
    "    nGenCarac = pd.DataFrame([genCarac['Characteristics ID'].iloc[idx] for idx,_ in enumerate(genCarac['Characteristics ID']) if idx+1 not in genCarac['Characteristics ID'] or type(genCarac['Component'].iloc[idx+1]) is int],index=genCarac['Component'].dropna().unique(),columns=['Number Caracteristics'])\n",
    "    Genlist = [Gen(0,0,i,0,0,0,0,0) for i in nGenCarac.index]\n",
    "    aux_nGenCarac = pd.DataFrame(np.cumsum(nGenCarac.values),index=nGenCarac.index)\n",
    "        \n",
    "    for i in nGenCarac.index: \n",
    "        if i-1 not in nGenCarac.index:\n",
    "            j=0\n",
    "            jlim = (aux_nGenCarac.loc[i].values-1).item()\n",
    "            while j<=jlim:\n",
    "                if genCarac['Characteristics D.'].iloc[j]=='Bus Location':\n",
    "                    Genlist[int(i-1)].bus = genCarac['Characteristics V.'].iloc[j]\n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Droop Contribution':\n",
    "                    Genlist[int(i-1)].dp_cont = genCarac['Characteristics V.'].iloc[j]\n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Gen Type':                    \n",
    "                    Genlist[i-1].gentype = genCarac['Characteristics V.'].iloc[j]                     \n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Nominal Voltage (V)':\n",
    "                    Genlist[i-1].vn = genCarac['Characteristics V.'].iloc[j]\n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Genlist[i-1].pn = genCarac['Characteristics V.'].iloc[j]\n",
    "                j+=1\n",
    "        else:\n",
    "            j = (aux_nGenCarac.loc[i-1].values).item()\n",
    "            jlim = (aux_nGenCarac.loc[i].values-1).item()\n",
    "            while j<=jlim:\n",
    "                if genCarac['Characteristics D.'].iloc[j]=='Bus Location':\n",
    "                    Genlist[i-1].bus = genCarac['Characteristics V.'].iloc[j]\n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Droop Contribution':\n",
    "                    Genlist[i-1].dp_cont = genCarac['Characteristics V.'].iloc[j]\n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Gen Type':                    \n",
    "                    Genlist[i-1].gentype = genCarac['Characteristics V.'].iloc[j]                    \n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Nominal Voltage (V)':\n",
    "                    Genlist[i-1].vn = genCarac['Characteristics V.'].iloc[j]\n",
    "                elif genCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Genlist[i-1].pn = genCarac['Characteristics V.'].iloc[j]\n",
    "                j+=1\n",
    "\n",
    "    #print(genCarac)\n",
    "\n",
    "    # Get a matrix of all voltages and powers points of Current/OS voltage/power droop curves from gen_data DataFrame to perform the interpolations\n",
    "\n",
    "    aux_namCol = ['Gen ID','Unnamed: 4','Total Time (h)']\n",
    "    for i,j in enumerate(gen_data.columns):\n",
    "        if type(j) is int:                       \n",
    "            aux_namCol.append(j)                    # Gives the correct column names to search the droop information in time in gen_data DataFrame \n",
    "\n",
    "\n",
    "    genDroop = gen_data[aux_namCol] \n",
    "    genDroop = genDroop.rename(columns={'Gen ID':genDroop['Gen ID'][1],'Unnamed: 4':genDroop['Unnamed: 4'][1],'Total Time (h)':genDroop['Total Time (h)'][1]})\n",
    "    genDroop = genDroop[2:]\n",
    "    maxDroopsGen = max(genDroop['Droop ID'].unique())\n",
    "    genDroop.index = range(0,genDroop.shape[0])\n",
    "\n",
    "    #print(genDroop)\n",
    "\n",
    "    for i in nGenCarac.index:\n",
    "        aux_indx = (np.array(genDroop.index[genDroop['Component']==i].tolist())).item()\n",
    "        if aux_indx is not None:\n",
    "            for j in range(0,len(Genlist)):\n",
    "                if Genlist[j].gen_id == i:\n",
    "                    Genlist[j].droop = genDroop.loc[aux_indx:aux_indx+maxDroopsGen-1]\n",
    "\n",
    "    for i in range(0,len(Genlist)):                # Sorting process (generator by generator)\n",
    "        \n",
    "        auxV = pd.DataFrame(Genlist[i].droop.loc[Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains('V')],'Droop'],index=Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains('V')].tolist())      # List of the names of the voltage points of the Current/OS voltage droop control curves of each generator\n",
    "        auxV.sort_values(by='Droop')                                                                                                                                                                                        # This command guarrantees that the voltage names are sorted\n",
    "\n",
    "        auxP = pd.DataFrame(Genlist[i].droop.loc[Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)],'Droop'],index=Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)].tolist())      # List of the names of the voltage points of the Current/OS voltage droop control curves of each generator\n",
    "        auxP.sort_values(by='Droop')                                                                                                                                                                                        # This command guarrantees that the voltage names are sorted\n",
    "\n",
    "        auxVP = list(zip(auxV.index,auxP.index))    # List with tuples with the indexes corresponding to the names of the \"coordinates\" of the points of the Current/OS voltage droop control curves of each generator\n",
    "\n",
    "        for k in Genlist[i].droop.columns:\n",
    "            if type(k) is int:\n",
    "                droop_list = []\n",
    "                for j,_ in enumerate(auxVP): \n",
    "                    droop_list.append(Genlist[i].droop.loc[np.array(list(auxVP[j])),k].values)  \n",
    "                df = pd.DataFrame(droop_list, columns=['x','y'])\n",
    "                df_sorted = df.sort_values(by=['x','y'])\n",
    "                Genlist[i].droop.loc[auxV.index,k] = df_sorted['x'].values\n",
    "                Genlist[i].droop.loc[auxP.index,k] = df_sorted['y'].values\n",
    "\n",
    "\n",
    "    #print(Genlist)\n",
    "\n",
    "    # Storage equipment \n",
    "    #print(\">Display Storage data\") \n",
    "    #print(stor_data)\n",
    "\n",
    "    # Get the caracteristics of the storage components in grid from stor_data DataFrame\n",
    "\n",
    "    storCarac = stor_data[['Storage ID','Unnamed: 1','Unnamed: 2','Unnamed: 3']]\n",
    "    storCarac = storCarac.rename(columns={'Storage ID':storCarac['Storage ID'][1],'Unnamed: 1':storCarac['Unnamed: 1'][1],'Unnamed: 2':storCarac['Unnamed: 2'][1],'Unnamed: 3':storCarac['Unnamed: 3'][1]})\n",
    "    storCarac = storCarac[2:]\n",
    "    storCarac.index = range(0,len(storCarac['Component']))\n",
    "    aux_idx = []\n",
    "    nStor = max(storCarac['Component'].unique())\n",
    "    for i,_ in enumerate(storCarac.iterrows()):\n",
    "        if (pd.isna(storCarac.iloc[i]).all()):\n",
    "            aux_idx.append(i)\n",
    "    storCarac = storCarac.drop(aux_idx)\n",
    "    storCarac.index = range(0,storCarac.shape[0])\n",
    "    nStorCarac = pd.DataFrame([storCarac['Characteristics ID'].iloc[idx] for idx,_ in enumerate(storCarac['Characteristics ID']) if idx+1 not in storCarac['Characteristics ID'] or type(storCarac['Component'].iloc[idx+1]) is int],index=storCarac['Component'].dropna().unique(),columns=['Number Caracteristics'])\n",
    "    Storlist = [Stor(0, 0, i, 0, 0, 0, 0, 0, 0, 0, 0, 0.00001, 0, 0) for i in nStorCarac.index]\n",
    "    aux_nStorCarac = pd.DataFrame(np.cumsum(nStorCarac.values),index=nStorCarac.index)\n",
    "        \n",
    "    for i in nStorCarac.index: \n",
    "        if i-1 not in nStorCarac.index:\n",
    "            j=0\n",
    "            jlim = aux_nStorCarac.loc[i].values-1\n",
    "            while j<=jlim:\n",
    "                if storCarac['Characteristics D.'].iloc[j]=='Bus Location':\n",
    "                    Storlist[int(i-1)].bus = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Droop Contribution':\n",
    "                    Storlist[i-1].dp_cont = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Initial SOC (%)':\n",
    "                    Storlist[int(i-1)].ini_SOC = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Storage Type':                    \n",
    "                    Storlist[i-1].gentype = storCarac['Characteristics V.'].iloc[j]                     \n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Nominal Voltage (V)':\n",
    "                    Storlist[i-1].vn = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Storlist[i-1].pn = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Maximum Energy (kWh)':\n",
    "                    Storlist[i-1].max_ener = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Maximum Power Delivery (kW)':\n",
    "                    Storlist[i-1].max_p_discharge = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Maximum Power Charge (kW)':\n",
    "                    Storlist[i-1].max_p_charge = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Upper SOC limit (%)':\n",
    "                    Storlist[i-1].max_SOC = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Lower SOC limit (%)':\n",
    "                    Storlist[i-1].min_SOC = storCarac['Characteristics V.'].iloc[j]\n",
    "                j+=1\n",
    "        else:\n",
    "            j=int(aux_nStorCarac.loc[i-1].values)\n",
    "            jlim = int(aux_nStorCarac.loc[i].values-1)\n",
    "            while j<=jlim:\n",
    "                if storCarac['Characteristics D.'].iloc[j]=='Bus Location':\n",
    "                    Storlist[int(i-1)].bus = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Droop Contribution':\n",
    "                    Storlist[i-1].dp_cont = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Initial SOC (%)':\n",
    "                    Storlist[int(i-1)].ini_SOC = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Storage Type':                    \n",
    "                    Storlist[i-1].gentype = storCarac['Characteristics V.'].iloc[j]                     \n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Nominal Voltage (V)':\n",
    "                    Storlist[i-1].vn = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Nominal Power (kW)':\n",
    "                    Storlist[i-1].pn = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Maximum Energy (kWh)':\n",
    "                    Storlist[i-1].max_ener = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Maximum Power Delivery (kW)':\n",
    "                    Storlist[i-1].max_p_discharge = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Maximum Power Charge (kW)':\n",
    "                    Storlist[i-1].max_p_charge = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Upper SOC limit (%)':\n",
    "                    Storlist[i-1].max_SOC = storCarac['Characteristics V.'].iloc[j]\n",
    "                elif storCarac['Characteristics D.'].iloc[j]=='Lower SOC limit (%)':\n",
    "                    Storlist[i-1].min_SOC = storCarac['Characteristics V.'].iloc[j]\n",
    "                j+=1\n",
    "\n",
    "    #print(storCarac)\n",
    "\n",
    "    # Get a matrix of all voltages and powers points of Current/OS voltage/power droop curves from stor_data DataFrame to perform the interpolations\n",
    "\n",
    "        # Get a matrix of all voltages and powers points of Current/OS voltage/power droop curves from gen_data DataFrame to perform the interpolations\n",
    "\n",
    "    aux_namCol = ['Storage ID','Unnamed: 4','Total Time (h)']\n",
    "    for i,j in enumerate(stor_data.columns):\n",
    "        if type(j) is int:                       \n",
    "            aux_namCol.append(j)                    # Gives the correct column names to search the droop information in time in gen_data DataFrame \n",
    "\n",
    "\n",
    "    storDroop = stor_data[aux_namCol] \n",
    "    storDroop = storDroop.rename(columns={'Storage ID':storDroop['Storage ID'][1],'Unnamed: 4':storDroop['Unnamed: 4'][1],'Total Time (h)':storDroop['Total Time (h)'][1]})\n",
    "    storDroop = storDroop[2:]\n",
    "    maxDroopsStor = max(storDroop['Droop ID'].unique())\n",
    "    storDroop.index = range(0,storDroop.shape[0])\n",
    "\n",
    "    #print(genDroop)\n",
    "\n",
    "    for i in nStorCarac.index:\n",
    "        aux_indx = (np.array(storDroop.index[storDroop['Component']==i].tolist())).item()\n",
    "        if aux_indx is not None:\n",
    "            for j in range(0,len(Storlist)):\n",
    "                if Storlist[j].stor_id == i:\n",
    "                    Storlist[j].droop = storDroop.loc[aux_indx:aux_indx+maxDroopsStor-1]\n",
    "\n",
    "    for i in range(0,len(Storlist)):                # Sorting process (generator by generator)\n",
    "        \n",
    "        auxV = pd.DataFrame(Storlist[i].droop.loc[Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains('V')],'Droop'],index=Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains('V')].tolist())      # List of the names of the voltage points of the Current/OS voltage droop control curves of each generator\n",
    "        auxV.sort_values(by='Droop')                                                                                                                                                                                        # This command guarrantees that the voltage names are sorted\n",
    "\n",
    "        auxP = pd.DataFrame(Storlist[i].droop.loc[Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)],'Droop'],index=Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)].tolist())      # List of the names of the voltage points of the Current/OS voltage droop control curves of each generator\n",
    "        auxP.sort_values(by='Droop')                                                                                                                                                                                        # This command guarrantees that the voltage names are sorted\n",
    "\n",
    "        auxVP = list(zip(auxV.index,auxP.index))    # List with tuples with the indexes corresponding to the names of the \"coordinates\" of the points of the Current/OS voltage droop control curves of each generator\n",
    "\n",
    "        for k in Storlist[i].droop.columns:\n",
    "            if type(k) is int:\n",
    "                droop_list = []\n",
    "                for j,_ in enumerate(auxVP): \n",
    "                    droop_list.append(Storlist[i].droop.loc[np.array(list(auxVP[j])),k].values)  \n",
    "                df = pd.DataFrame(droop_list, columns=['x','y'])\n",
    "                df_sorted = df.sort_values(by=['x','y'])\n",
    "                Storlist[i].droop.loc[auxV.index,k] = df_sorted['x'].values\n",
    "                Storlist[i].droop.loc[auxP.index,k] = df_sorted['y'].values\n",
    "\n",
    "\n",
    "    #print(Storlist)\n",
    "\n",
    "    ########## NETWORK CREATION ###########\n",
    "\n",
    "    # Create an empty pandapower network\n",
    "    #print(\">Create an empty pandapower network\")\n",
    "    net = pp.create_empty_network(name='', f_hz=50.0, sn_mva=base_power.item()/1e3)\n",
    "    #print(net)\n",
    "\n",
    "    \n",
    "    # Create buses\n",
    "    #print(\">Create buses\") \n",
    "    \n",
    "    buses = {}\n",
    "    \n",
    "    for i in range(0,len(Nodelist)):\n",
    "        buses[Nodelist[i].node_number.item()] = pp.create_bus(net, vn_kv=Nodelist[i].vn_kv.item(), name=f\"Bus {Nodelist[i].node_name.item()}\") # Convert from V to kV (Important to PF convergence)\n",
    "    \n",
    "    #print(buses)\n",
    "\n",
    "    for i in range(0,len(Linelist)):\n",
    "\n",
    "        pp.create_line_from_parameters(net,\n",
    "                                    name='Line ' + str(Linelist[i].line_id),\n",
    "                                    from_bus=buses[Linelist[i].from_bus.item()],\n",
    "                                    to_bus=buses[Linelist[i].to_bus.item()],\n",
    "                                    length_km=Linelist[i].l_km.item(),\n",
    "                                    r_ohm_per_km=Linelist[i].r_ohm_km.item(),\n",
    "                                    x_ohm_per_km=Linelist[i].x_ohm_km.item(),\n",
    "                                    c_nf_per_km=Linelist[i].c_ohm_km.item(),\n",
    "                                    max_i_ka=Linelist[i].Imax_kA.item())\n",
    "\n",
    "    # Create transformers\n",
    "    #print(\">Create transformers\") \n",
    "\n",
    "    for i in range(0,len(Trafolist)):\n",
    "\n",
    "        aux_fbus_idx = 0\n",
    "        aux_tbus_idx = 0\n",
    "\n",
    "        for j in range(0,len(Nodelist)):\n",
    "            if Nodelist[j].node_number.item() == Trafolist[i].from_bus.item():\n",
    "                aux_fbus_idx = j\n",
    "            if Nodelist[j].node_number.item() == Trafolist[i].to_bus.item():\n",
    "                aux_tbus_idx = j\n",
    "\n",
    "        pp.create_std_type(net, {\"sn_mva\": Trafolist[i].sn.item(),\n",
    "            \"vn_hv_kv\": Nodelist[aux_fbus_idx].vn_kv.item(),\n",
    "            \"vn_lv_kv\": Nodelist[aux_tbus_idx].vn_kv.item(),\n",
    "            \"vk_percent\": Trafolist[i].vk_percent.item(),\n",
    "            \"vkr_percent\": Trafolist[i].vkr_percent.item(),\n",
    "            \"pfe_kw\": Trafolist[i].pfe_kw.item(),\n",
    "            \"i0_percent\": Trafolist[i].i0_percent.item(),\n",
    "            \"shift_degree\": Trafolist[i].shift_degree.item(),\n",
    "            \"vector_group\": str(Trafolist[i].vector_group.item()),\n",
    "            \"tap_side\": str(Trafolist[i].tap_side.item()),\n",
    "            \"tap_neutral\": 0, # Miss in Excel file\n",
    "            \"tap_min\": Trafolist[i].min_tap.item(),\n",
    "            \"tap_max\": Trafolist[i].max_tap.item(),\n",
    "            \"tap_step_degree\": Trafolist[i].tap_step_degree.item(),\n",
    "            \"tap_step_percent\": Trafolist[i].tap_step.item(),\n",
    "            \"tap_phase_shifter\": Trafolist[i].tap_ps.item()=='True',\n",
    "            \"vk0_percent\": Trafolist[i].vk0_percent.item(),\n",
    "            \"vkr0_percent\": Trafolist[i].vkr0_percent,\n",
    "            \"mag0_percent\": Trafolist[i].mag0_percent.item(),\n",
    "            \"mag0_rx\": Trafolist[i].mag0_rx.item(),\n",
    "            \"si0_hv_partial\": Trafolist[i].si0_hv_partial.item(),}, name='Transformer', element=\"trafo\")\n",
    "\n",
    "        pp.create_transformer(net, hv_bus = buses[Trafolist[i].from_bus.item()], lv_bus = buses[Trafolist[i].to_bus.item()], name = 'Transformer ' + str(Trafolist[i].trafo_id), std_type = 'Transformer')\n",
    "\n",
    "    # Create converters\n",
    "    #print(\">Create converters\") \n",
    "\n",
    "    for i in range(0,len(Convlist)):\n",
    "\n",
    "        if Convlist[i].convtype == 'Converter DC/DC':\n",
    "\n",
    "            # For the DC/DC converters we assume a analogy with three phase transformers (pandapower model -> transformer model)\n",
    "\n",
    "            aux_fbus_idx = 0\n",
    "            aux_tbus_idx = 0\n",
    "\n",
    "            for j in range(0,len(Nodelist)):\n",
    "                if Nodelist[j].node_number.item() == Convlist[i].from_bus:\n",
    "                    aux_fbus_idx = j\n",
    "                if Nodelist[j].node_number.item() == Convlist[i].to_bus:\n",
    "                    aux_tbus_idx = j\n",
    "\n",
    "            pp.create_std_type(net, {\"sn_mva\": Convlist[i].pn,\n",
    "                \"vn_hv_kv\": Nodelist[aux_fbus_idx].vn_kv.item(),\n",
    "                \"vn_lv_kv\": Nodelist[aux_tbus_idx].vn_kv.item(),\n",
    "                \"vk_percent\": 0.01,\n",
    "                \"vkr_percent\": 0.01,\n",
    "                \"pfe_kw\": 0,\n",
    "                \"i0_percent\": 0,\n",
    "                \"shift_degree\": 0,\n",
    "                \"vector_group\": 'Yyn',\n",
    "                \"tap_side\": 'lv',\n",
    "                \"tap_neutral\": 0, \n",
    "                \"tap_min\":-2,\n",
    "                \"tap_max\": 2,\n",
    "                \"tap_step_degree\": 0,\n",
    "                \"tap_step_percent\": 2.5,\n",
    "                \"tap_phase_shifter\": False,\n",
    "                \"vk0_percent\": 0,\n",
    "                \"vkr0_percent\": 0,\n",
    "                \"mag0_percent\": 0,\n",
    "                \"mag0_rx\": 0,\n",
    "                \"si0_hv_partial\": 0,}, name='Converter DC/DC', element=\"trafo\")\n",
    "\n",
    "            pp.create_transformer(net, hv_bus = buses[Convlist[i].from_bus], lv_bus = buses[Convlist[i].to_bus], name = 'Converter DC/DC ' + str(Convlist[i].conv_id), std_type = 'Converter DC/DC')\n",
    "\n",
    "        elif Convlist[i].convtype == 'Converter AC/DC': \n",
    "\n",
    "            # Also for the AC/DC converters we assume a two network element model (pandapower models -> transformer model + ward model)\n",
    "\n",
    "            aux_fbus_idx = 0\n",
    "            aux_tbus_idx = 0\n",
    "\n",
    "            for j in range(0,len(Nodelist)):\n",
    "                \n",
    "                # Transformer buses automatic finding\n",
    "                \n",
    "                if Nodelist[j].node_number.item() == Convlist[i].from_bus:\n",
    "                    aux_fbus_idx = j\n",
    "                if Nodelist[j].node_number.item() == Convlist[i].to_bus:\n",
    "                    aux_tbus_idx = j\n",
    "                \n",
    "                # Ward bus automatic finding\n",
    "                \n",
    "                if Nodelist[j].elec_type.item() == 'AC':\n",
    "                    aux_ward_idx = j\n",
    "                    if Nodelist[j].node_number.item() == Convlist[i].from_bus:\n",
    "                        aux_ward_bus_idx = Convlist[i].from_bus\n",
    "                    elif Nodelist[j].node_number.item() == Convlist[i].to_bus:\n",
    "                        aux_ward_bus_idx = Convlist[i].to_bus\n",
    "\n",
    "\n",
    "            pp.create_std_type(net, {\"sn_mva\": Convlist[i].pn,\n",
    "                \"vn_hv_kv\": Nodelist[aux_fbus_idx].vn_kv.item(),\n",
    "                \"vn_lv_kv\": Nodelist[aux_tbus_idx].vn_kv.item(),\n",
    "                \"vk_percent\": 0.01,\n",
    "                \"vkr_percent\": 0.009,\n",
    "                \"pfe_kw\": 0,\n",
    "                \"i0_percent\": 0,\n",
    "                \"shift_degree\": 0,\n",
    "                \"vector_group\": 'Yyn',\n",
    "                \"tap_side\": 'lv',\n",
    "                \"tap_neutral\": 0, \n",
    "                \"tap_min\":-2,\n",
    "                \"tap_max\": 2,\n",
    "                \"tap_step_degree\": 0,\n",
    "                \"tap_step_percent\": 2.5,\n",
    "                \"tap_phase_shifter\": False,\n",
    "                \"vk0_percent\": 0,\n",
    "                \"vkr0_percent\": 0,\n",
    "                \"mag0_percent\": 0,\n",
    "                \"mag0_rx\": 0,\n",
    "                \"si0_hv_partial\": 0,}, name='Converter AC/DC', element=\"trafo\")\n",
    "\n",
    "            pp.create_transformer(net, hv_bus = buses[Convlist[i].from_bus], lv_bus = buses[Convlist[i].to_bus], name = 'Converter AC/DC ' + str(Convlist[i].conv_id), std_type = 'Converter AC/DC')    # Creation of transformer element\n",
    "\n",
    "            pp.create_ward(net, buses[aux_ward_bus_idx], 0, 0, 0, 0, name = 'Converter AC/DC ' + str(Convlist[i].conv_id))  # Creation of ward element and corresponding initialization of reactive power in zero\n",
    "\n",
    "            #pp.create_vsc(net,                              # It is dependent from vsc development from pandapower\n",
    "            #    bus = buses[Convlist[i].from_bus],\n",
    "            #    bus_dc = buses[Convlist[i].to_bus],\n",
    "            #    r_ohm = Convlist[i].r_ohm.item(),\n",
    "            #    x_ohm = Convlist[i].x_ohm.item(),\n",
    "            #    r_dc_ohm = Convlist[i].r_dc_ohm.item(),\n",
    "            #    pl_dc_mw = Convlist[i].pl_dc_mw.item(),\n",
    "            #    control_mode_ac = str(Convlist[i].control_mode_ac.item()),  \n",
    "            #    control_value_ac = Convlist[i].control_value_ac.item(), \n",
    "            #    control_mode_dc = str(Convlist[i].control_mode_dc.item()), \n",
    "            #    control_value_dc = Convlist[i].control_value_dc.item(),\n",
    "            #    name='Converter AC/DC ' + str(Convlist[i].conv_id))\n",
    "\n",
    "    # Create loads, gens and storage (in progress)\n",
    "    #print(\">Create loads, generators and storage\")\n",
    "\n",
    "    # Load equipment\n",
    "    #print(\">Create loads\")\n",
    "\n",
    "    for i in range(0,len(Loadlist)):\n",
    "        auxV = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains('V')]\n",
    "        auxP = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)]\n",
    "        auxPS = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'P[A-Za-z]+', regex=True)]\n",
    "        auxQS = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'Q[A-Za-z]+', regex=True)]\n",
    "        auxT = min([i for i in Loadlist[i].droop.columns if type(i) is int])\n",
    "        for j in range(0,len(Nodelist)):\n",
    "            if Nodelist[j].node_number.item() == Loadlist[i].bus:\n",
    "                aux_electype = Nodelist[j].elec_type.item()\n",
    "        if aux_electype == 'DC':\n",
    "            if Loadlist[i].dp_cont == 1:\n",
    "                pdroop = np.interp(1,Loadlist[i].droop.loc[auxV,auxT].values,Loadlist[i].droop.loc[auxP,auxT].values)\n",
    "                pset = Loadlist[i].droop.loc[auxPS,auxT].values\n",
    "                pload = 0\n",
    "                if abs(pset) > abs(pdroop): # Verification if setpoint of power is higher than droop curve\n",
    "                    pload = pdroop\n",
    "                else:\n",
    "                    pload = pset\n",
    "\n",
    "                pp.create_load(net,\n",
    "                            name=Loadlist[i].load_id,\n",
    "                            bus=buses[Loadlist[i].bus], # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=pload*Loadlist[i].pn/1000,  # Convert kW to MW \n",
    "                            q_mvar=0)  # assuming zero reactive power\n",
    "            else:\n",
    "                pp.create_load(net,\n",
    "                            name=Loadlist[i].load_id,\n",
    "                            bus=buses[Loadlist[i].bus], # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=(Loadlist[i].droop.loc[auxPS,auxT])*Loadlist[i].pn/1000,  # Convert kW to MW, it is assumed the maximum power point as constant power load situation (dp_cont variable equal to 0)\n",
    "                            q_mvar=0)  # assuming zero reactive power\n",
    "        elif aux_electype == 'AC':\n",
    "            pp.create_load(net,\n",
    "                            name=Loadlist[i].load_id,\n",
    "                            bus=buses[Loadlist[i].bus], # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=(Loadlist[i].droop.loc[auxPS,auxT])*Loadlist[i].pn/1000,  # Convert kW to MW, it is assumed the maximum power point as constant power load situation (dp_cont variable equal to 0)\n",
    "                            q_mvar=(Loadlist[i].droop.loc[auxQS,auxT])*Loadlist[i].pn/1000)  # assuming zero reactive power\n",
    "        \n",
    "    # Generator equipment\n",
    "    #print(\">Create generators\")\n",
    "\n",
    "    for i in range(0,len(Genlist)):\n",
    "        auxV = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains('V')]\n",
    "        auxP = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)]\n",
    "        auxPS = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'P[A-Za-z]+', regex=True)]\n",
    "        auxQS = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'Q[A-Za-z]+', regex=True)]\n",
    "        auxT = min([i for i in Genlist[i].droop.columns if type(i) is int])\n",
    "        for j in range(0,len(Nodelist)):\n",
    "            if Nodelist[j].node_number.item() == Genlist[i].bus:\n",
    "                aux_electype = Nodelist[j].elec_type.item()\n",
    "        if aux_electype == 'DC':\n",
    "            if Genlist[i].gen_id == slack_comp.values:\n",
    "                if Genlist[i].dp_cont == 1:\n",
    "                    pdroop = np.interp(1,Genlist[i].droop.loc[auxV,auxT].values,Genlist[i].droop.loc[auxP,auxT].values)\n",
    "                    pset = Genlist[i].droop.loc[auxPS,auxT].values\n",
    "                    pload = 0\n",
    "\n",
    "                    if abs(pset) > abs(pdroop): # Verification if setpoint of power is higher than droop curve\n",
    "                        pload = pdroop\n",
    "                    else:\n",
    "                        pload = pset\n",
    "\n",
    "                    pp.create_gen(net,\n",
    "                            name=Genlist[i].gen_id,\n",
    "                            bus=buses[Genlist[i].bus],      # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=pset*Genlist[i].pn/1000,   # Convert kW to MW \n",
    "                            vm_pu=1,                        # assuming a voltage setpoint of 1 (Can be an improvement to make by saying that one of the caracteristics of the gen is the voltage setpoint)\n",
    "                            slack=True)  \n",
    "                else:\n",
    "                    pp.create_gen(net,\n",
    "                            name=Genlist[i].gen_id,\n",
    "                            bus=buses[Genlist[i].bus], \n",
    "                            p_mw=(Genlist[i].droop.loc[auxPS,auxT])*Genlist[i].pn/1000,  # Convert kW to MW, it is assumed the maximum power point as constant power gen situation (dp_cont variable equal to 0)\n",
    "                            vm_pu=1,\n",
    "                            slack=True)      \n",
    "            elif Genlist[i].dp_cont == 1: # Non-slack generators (These generators must be static ones since it is only known the P or by its output power or by its droop curve)\n",
    "                pdroop = np.interp(1,Genlist[i].droop.loc[auxV,auxT].values,Genlist[i].droop.loc[auxP,auxT].values)\n",
    "                pset = Genlist[i].droop.loc[auxPS,auxT].values\n",
    "                pload = 0\n",
    "                if abs(pset) > abs(pdroop): # Verification if setpoint of power is higher than droop curve\n",
    "                    pload = pdroop\n",
    "                else:\n",
    "                    pload = pset\n",
    "                pp.create_sgen(net,\n",
    "                            name=Genlist[i].gen_id,\n",
    "                            bus=buses[Genlist[i].bus], \n",
    "                            p_mw=np.interp(1,Genlist[i].droop.loc[auxV,auxT].values,Genlist[i].droop.loc[auxP,auxT].values)*Genlist[i].pn/1000)  \n",
    "            else:\n",
    "                pp.create_sgen(net,\n",
    "                            name=Genlist[i].gen_id,\n",
    "                            bus=buses[Genlist[i].bus], \n",
    "                            p_mw=(Genlist[i].droop.loc[auxPS,auxT])*Genlist[i].pn/1000)   \n",
    "        elif aux_electype == 'AC':\n",
    "            if Genlist[i].gen_id == slack_comp.values:\n",
    "                pp.create_gen(net,\n",
    "                            name=Genlist[i].gen_id,\n",
    "                            bus=buses[Genlist[i].bus], # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=(Genlist[i].droop.loc[auxPS,auxT])*Genlist[i].pn/1000,  # Convert kW to MW, it is assumed the maximum power point as constant power load situation (dp_cont variable equal to 0)\n",
    "                            vm_pu=1,\n",
    "                            max_q_mvar=(Genlist[i].droop.loc[auxQS,auxT])*Genlist[i].pn/1000,\n",
    "                            slack = True)  \n",
    "            else:\n",
    "                pp.create_sgen(net,\n",
    "                            name=Genlist[i].gen_id,\n",
    "                            bus=buses[Genlist[i].bus], \n",
    "                            p_mw=(Genlist[i].droop.loc[auxPS,auxT])*Genlist[i].pn/1000,\n",
    "                            max_q_mvar=(Genlist[i].droop.loc[auxQS,auxT])*Genlist[i].pn/1000)   \n",
    "\n",
    "    # Storage equipment\n",
    "    #print(\">Create BESS\")  # Misses how to consider a BESS element to be the slack component\n",
    "\n",
    "    for i in range(0,len(Storlist)):\n",
    "        auxV = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains('V')]\n",
    "        auxP = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)]\n",
    "        auxPS = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'P[A-Za-z]+', regex=True)]\n",
    "        auxQS = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'Q[A-Za-z]+', regex=True)]\n",
    "        auxT = min([i for i in Storlist[i].droop.columns if type(i) is int])\n",
    "        for j in range(0,len(Nodelist)):\n",
    "            if Nodelist[j].node_number.item() == Storlist[i].bus:\n",
    "                aux_electype = Nodelist[j].elec_type.item()\n",
    "        if aux_electype == 'DC':\n",
    "            if Storlist[i].dp_cont == 1:\n",
    "                pdroop = np.interp(1,Storlist[i].droop.loc[auxV,auxT].values,Storlist[i].droop.loc[auxP,auxT].values)\n",
    "                pset = Storlist[i].droop.loc[auxPS,auxT].values\n",
    "                pload = 0\n",
    "                if abs(pset) > abs(pdroop): # Verification if setpoint of power is higher than droop curve\n",
    "                    pload = pdroop\n",
    "                else:\n",
    "                    pload = pset\n",
    "\n",
    "                pp.create_storage(net,\n",
    "                            name=Storlist[i].stor_id,\n",
    "                            bus=buses[Storlist[i].bus], # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=pload*Storlist[i].pn/1000,  # Convert kW to MW \n",
    "                            q_mvar=0,    # assuming zero reactive power (DC)\n",
    "                            max_e_mwh=Storlist[i].max_ener/1000, # Convert from kWh to MWh\n",
    "                            soc_percent=Storlist[i].ini_SOC)   \n",
    "            else:\n",
    "                pp.create_storage(net,\n",
    "                            name=Storlist[i].stor_id,\n",
    "                            bus=buses[Storlist[i].bus], # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=(Storlist[i].droop.loc[auxPS,auxT])*Storlist[i].pn/1000,  # Convert kW to MW, it is assumed the maximum power point as constant power load situation (dp_cont variable equal to 0)\n",
    "                            q_mvar=0,\n",
    "                            max_e_mwh=Storlist[i].max_ener/1000, # Convert from kWh to MWh\n",
    "                            soc_percent=Storlist[i].ini_SOC)\n",
    "        elif aux_electype == 'AC':\n",
    "            pp.create_storage(net,\n",
    "                            name=Storlist[i].stor_id,\n",
    "                            bus=buses[Storlist[i].bus], # The interpolation is the initial droop adjustment assuming a flat start in PF solving\n",
    "                            p_mw=(Storlist[i].droop.loc[auxPS,auxT])*Storlist[i].pn/1000,  # Convert kW to MW, it is assumed the maximum power point as constant power load situation (dp_cont variable equal to 0)\n",
    "                            q_mvar=(Storlist[i].droop.loc[auxQS,auxT])*Storlist[i].pn/1000,  # assuming zero reactive power\n",
    "                            max_e_mwh=Storlist[i].max_ener/1000, # Convert from kWh to MWh\n",
    "                            soc_percent=Storlist[i].ini_SOC)\n",
    "\n",
    "    return buses, base_power, Nodelist, Linelist, Trafolist, Convlist, Loadlist, Genlist, Storlist, net, period_duration, slack_comp\n",
    "\n",
    "[buses, base_power, Nodelist, Linelist, Trafolist, Convlist, Loadlist, Genlist, Storlist, net, period_duration, slack_comp] = Griddata(general_data,node_data,line_data,transformer_data,converter_data,load_data,gen_data,stor_data)\n",
    "\n",
    "\n",
    "# Misses how to handle with empty lists, for exemple, in case that are not used transformers in network\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plotting the droop curves of each load and generator equipment \n",
    "\n",
    "# Load equipment\n",
    "# print('\\nPlotting the VP droop curve caracteristic of each load')\n",
    "\n",
    "print(\">Plot droop evolution\")\n",
    "\n",
    "time_list = [i for i in Loadlist[0].droop.columns if type(i) is int]\n",
    "\n",
    "fig, axs = plt.subplots(len(time_list),1,figsize=(14, 30),layout='constrained')\n",
    "#fig.title(\"Droop Evolution\")\n",
    "labels = []\n",
    "for i in range(0,len(time_list)):\n",
    "    for j in range(0,len(Loadlist)):\n",
    "        auxV = Loadlist[j].droop.index[Loadlist[j].droop['Droop'].str.contains('V')]\n",
    "        auxP = Loadlist[j].droop.index[Loadlist[j].droop['Droop'].str.contains(r'P\\d+', regex=True)]   \n",
    "        axs[i].plot(Loadlist[j].droop.loc[auxP,time_list[i]].values,Loadlist[j].droop.loc[auxV,time_list[i]].values, marker='o', linestyle='-')     \n",
    "        labels.append(['Load' + str(Loadlist[j].load_id)])\n",
    "    axs[i].set_ylabel(\"Voltage (pu)\")\n",
    "    axs[i].legend(labels)\n",
    "    axs[i].grid(True)\n",
    "    \n",
    "plt.show()\n",
    "\n",
    "# Generator equipment\n",
    "# print('\\nPlotting the VP droop curve caracteristic of each generator')\n",
    "\n",
    "print(\">Plot droop evolution\")\n",
    "\n",
    "fig, axs = plt.subplots(len(time_list),1,figsize=(14, 30),layout='constrained')\n",
    "#fig.title(\"Droop Evolution\")\n",
    "labels = []\n",
    "for i in range(0,len(time_list)):\n",
    "    for j in range(0,len(Genlist)):\n",
    "        auxV = Genlist[j].droop.index[Genlist[j].droop['Droop'].str.contains('V')]\n",
    "        auxP = Genlist[j].droop.index[Genlist[j].droop['Droop'].str.contains(r'P\\d+', regex=True)]   \n",
    "        axs[i].plot(Genlist[j].droop.loc[auxP,time_list[i]].values,Genlist[j].droop.loc[auxV,time_list[i]].values, marker='o', linestyle='-')     \n",
    "        labels.append(['Generator' + str(Genlist[j].gen_id)])\n",
    "    axs[i].set_ylabel(\"Voltage (pu)\")\n",
    "    axs[i].legend(labels)\n",
    "    axs[i].grid(True)\n",
    "    \n",
    "plt.show()\n",
    "\n",
    "\n",
    "\n",
    "# Storage equipment\n",
    "# print('\\nPlotting the VP droop curve caracteristic of each storage')\n",
    "\n",
    "print(\">Plot droop evolution\")\n",
    "\n",
    "fig, axs = plt.subplots(len(time_list),1,figsize=(14, 30),layout='constrained')\n",
    "#fig.title(\"Droop Evolution\")\n",
    "labels = []\n",
    "for i in range(0,len(time_list)):\n",
    "    for j in range(0,len(Storlist)):\n",
    "        auxV = Storlist[j].droop.index[Storlist[j].droop['Droop'].str.contains('V')]\n",
    "        auxP = Storlist[j].droop.index[Storlist[j].droop['Droop'].str.contains(r'P\\d+', regex=True)]   \n",
    "        axs[i].plot(Storlist[j].droop.loc[auxP,time_list[i]].values,Storlist[j].droop.loc[auxV,time_list[i]].values, marker='o', linestyle='-')     \n",
    "        labels.append(['Storage' + str(Storlist[j].stor_id)])\n",
    "    axs[i].set_ylabel(\"Voltage (pu)\")\n",
    "    axs[i].legend(labels)\n",
    "    axs[i].grid(True)\n",
    "    \n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Selection of Time-Variant or Steady-State PF\n",
    "\n",
    "def timeQuest(tl): # This function allows to give the pretended time series that the user want to study\n",
    "\n",
    "    tl_chosen = []\n",
    "\n",
    "    PF_type = input('\\nWhat type of Power Flow Analysis is pretended to make? For Time-Dependent, type TD, for Steady-State, type SS, to exit, EXIT')\n",
    "\n",
    "    if PF_type == 'TD':\n",
    "        tl_chosen = tl.copy()\n",
    "    elif PF_type == 'SS':\n",
    "        \n",
    "        print('\\nThe list of time instants that is present in the Excel is given by',tl) # Print the time data for user choose the time instant that want to study \n",
    "        usInput = input('\\nWhat is the time(s) instant(s) that the user want to study?') # Receive the hourly data from user (for exemple [1,2,...] or 1 2 ...) and treat that data\n",
    "        print(usInput)\n",
    "        if len(usInput) > 1:\n",
    "            aux_usInput = [tl[x] for x in range(0,len(tl)) if usInput in tl[x]]\n",
    "        else:\n",
    "            aux_usInput = usInput\n",
    "        tl_chosen.append(aux_usInput)\n",
    "    elif PF_type == 'EXIT':\n",
    "        return\n",
    "    else:\n",
    "        print('\\nCaution! The value introduced is not according to the expected input.')\n",
    "        return timeQuest(tl)\n",
    "\n",
    "    return tl_chosen\n",
    "\n",
    "def loaddroopReadjustment(net,Loadlist,bus_voltages,base_power,aux_Matrixp,T,t):\n",
    "\n",
    "    for i in range(0,len(Loadlist)): \n",
    "        \n",
    "        auxV = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains('V')]\n",
    "        auxP = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)]\n",
    "        auxPS = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'P[A-Za-z]+', regex=True)]\n",
    "        auxQS = Loadlist[i].droop.index[Loadlist[i].droop['Droop'].str.contains(r'Q[A-Za-z]+', regex=True)]\n",
    "        loadIdx = net.load.index[Loadlist[i].load_id==net.load['name']]\n",
    "\n",
    "        for j in range(0,len(Nodelist)):\n",
    "            if Nodelist[j].node_number.item() == Loadlist[i].bus:\n",
    "                aux_electype = Nodelist[j].elec_type.item()\n",
    "\n",
    "        if aux_electype == 'DC':\n",
    "\n",
    "            if Loadlist[i].dp_cont == 1:\n",
    "\n",
    "                pdroop = np.interp(bus_voltages[buses[Loadlist[i].bus]].item(),Loadlist[i].droop.loc[auxV,T].values,Loadlist[i].droop.loc[auxP,T].values)\n",
    "                pset = Loadlist[i].droop.loc[auxPS,T].values\n",
    "                pload = 0\n",
    "            \n",
    "                if abs(pset) > abs(pdroop): # Verification if setpoint of power is higher than droop curve\n",
    "                    pload = pdroop\n",
    "                else:\n",
    "                    pload = pset\n",
    "                \n",
    "                net.load.loc[loadIdx,'p_mw'] = pload*Loadlist[i].pn/1000\n",
    "            \n",
    "            elif Loadlist[i].dp_cont == 0 and t == 0:\n",
    "\n",
    "                net.load.loc[loadIdx,'p_mw']=(Loadlist[i].droop.loc[auxPS,T]).item()*Loadlist[i].pn/1000\n",
    "\n",
    "        elif aux_electype == 'AC':\n",
    "\n",
    "            net.load.loc[loadIdx,'p_mw']=(Loadlist[i].droop.loc[auxPS,T]).item()*Loadlist[i].pn/1000\n",
    "            net.load.loc[loadIdx,'q_mvar']=(Loadlist[i].droop.loc[auxQS,T]).item()*Loadlist[i].pn/1000\n",
    "\n",
    "        aux_Matrixp[i,t] = (float(net.load.loc[loadIdx,'p_mw'].iloc[0])*(1000/base_power.values)).item()\n",
    "\n",
    "    return net,aux_Matrixp\n",
    "\n",
    "def gendroopReadjustment(net,Genlist,bus_voltages,base_power,T):\n",
    "\n",
    "    for i in range(0,len(Genlist)): \n",
    "    \n",
    "        auxV = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains('V')]\n",
    "        auxP = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)]\n",
    "        auxPS = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'P[A-Za-z]+', regex=True)]\n",
    "        auxQS = Genlist[i].droop.index[Genlist[i].droop['Droop'].str.contains(r'Q[A-Za-z]+', regex=True)]\n",
    "        genIdx = net.gen.index[Genlist[i].gen_id==net.gen['name']] \n",
    "\n",
    "        for j in range(0,len(Nodelist)):\n",
    "            if Nodelist[j].node_number.item() == Genlist[i].bus:\n",
    "                aux_electype = Nodelist[j].elec_type.item()\n",
    "\n",
    "        if aux_electype == 'DC':\n",
    "\n",
    "            if Genlist[i].dp_cont == 1:\n",
    "\n",
    "                pdroop = np.interp(bus_voltages[buses[Genlist[i].bus]].item(),Genlist[i].droop.loc[auxV,T].values,Genlist[i].droop.loc[auxP,T].values)\n",
    "                pset = Genlist[i].droop.loc[auxPS,T].values\n",
    "                pgen = 0\n",
    "            \n",
    "                if abs(pset) > abs(pdroop): # Verification if setpoint of power is higher than droop curve\n",
    "                    pgen = pdroop\n",
    "                else:\n",
    "                    pgen = pset\n",
    "                \n",
    "                if Genlist[i].gen_id == slack_comp.values:\n",
    "                    net.gen.loc[genIdx,'p_mw'] = pgen*Genlist[i].pn/1000\n",
    "\n",
    "                else:\n",
    "                    net.sgen.loc[genIdx,'p_mw'] = pgen*Genlist[i].pn/1000\n",
    "            \n",
    "            elif Genlist[i].dp_cont == 0 and t == 0:\n",
    "\n",
    "                if Genlist[i].gen_id == slack_comp.values:\n",
    "\n",
    "                    net.gen.loc[genIdx,'p_mw'] = (Genlist[i].droop.loc[auxPS,T]).item()*Genlist[i].pn/1000\n",
    "\n",
    "                else:\n",
    "                    net.sgen.loc[genIdx,'p_mw'] = (Genlist[i].droop.loc[auxPS,T]).item()*Genlist[i].pn/1000\n",
    "\n",
    "            else:\n",
    "\n",
    "                if Genlist[i].gen_id == slack_comp.values:\n",
    "\n",
    "                    net.gen.loc[genIdx,'p_mw'] = (Genlist[i].droop.loc[auxPS,T]).item()*Genlist[i].pn/1000\n",
    "\n",
    "                else:\n",
    "                    net.sgen.loc[genIdx,'p_mw'] = (Genlist[i].droop.loc[auxPS,T]).item()*Genlist[i].pn/1000\n",
    "\n",
    "\n",
    "        elif aux_electype == 'AC':\n",
    "            \n",
    "            if Genlist[i].gen_id == slack_comp.values:\n",
    "                net.gen.loc[genIdx,'p_mw']=(Genlist[i].droop.loc[auxPS,T]).item()*Genlist[i].pn/1000\n",
    "                net.gen.loc[genIdx,'q_mvar']=(Genlist[i].droop.loc[auxQS,T]).item()*Genlist[i].pn/1000\n",
    "\n",
    "            else:\n",
    "                net.gen.loc[genIdx,'p_mw']=(Genlist[i].droop.loc[auxPS,T]).item()*Genlist[i].pn/1000\n",
    "                net.gen.loc[genIdx,'q_mvar']=(Genlist[i].droop.loc[auxQS,T]).item()*Genlist[i].pn/1000\n",
    "\n",
    "    return net\n",
    "\n",
    "def storControl(net,Storlist,bus_voltages,base_power,period_duration,aux_Matrixsoc,T,t): \n",
    "\n",
    "    for i in range(0,len(Storlist)):\n",
    "\n",
    "        auxV = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains('V')]\n",
    "        auxP = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'P\\d+', regex=True)]\n",
    "        auxPS = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'P[A-Za-z]+', regex=True)]\n",
    "        auxQS = Storlist[i].droop.index[Storlist[i].droop['Droop'].str.contains(r'Q[A-Za-z]+', regex=True)]\n",
    "        storIdx = net.storage.index[Storlist[i].stor_id==net.storage['name']]\n",
    "        P_stor = 0\n",
    "        Q_stor = 0\n",
    "\n",
    "        for j in range(0,len(Nodelist)):\n",
    "            if Nodelist[j].node_number.item() == Storlist[i].bus:\n",
    "                aux_electype = Nodelist[j].elec_type.item()\n",
    "\n",
    "        if aux_electype == 'DC':\n",
    "\n",
    "            if Storlist[i].dp_cont == 1:\n",
    "\n",
    "                pdroop = np.interp(bus_voltages[buses[Storlist[i].bus]].item(),Storlist[i].droop.loc[auxV,T].values,Storlist[i].droop.loc[auxP,T].values)\n",
    "                pset = Storlist[i].droop.loc[auxPS,T].values\n",
    "                pstor = 0\n",
    "            \n",
    "                if abs(pset) > abs(pdroop): # Verification if setpoint of power is higher than droop curve\n",
    "                    pstor = pdroop\n",
    "                else:\n",
    "                    pstor = pset\n",
    "                \n",
    "                net.storage.loc[storIdx,'p_mw'] = pstor*Storlist[i].pn/1000\n",
    "            \n",
    "            elif Storlist[i].dp_cont == 0 and t == 0:\n",
    "\n",
    "                net.storage.loc[storIdx,'p_mw']=(Storlist[i].droop.loc[auxPS,T])*Storlist[i].pn/1000\n",
    "            \n",
    "            P_stor = net.storage.loc[storIdx,'p_mw'].item()\n",
    "            Q_stor = 0\n",
    "\n",
    "        elif aux_electype == 'AC':\n",
    "\n",
    "            net.storage.loc[storIdx,'p_mw']=(Storlist[i].droop.loc[auxPS,T])*Storlist[i].pn/1000\n",
    "            net.storage.loc[storIdx,'q_mvar']=(Storlist[i].droop.loc[auxQS,T])*Storlist[i].pn/1000\n",
    "\n",
    "            P_stor = net.storage.loc[storIdx,'p_mw'].item()\n",
    "            Q_stor = net.storage.loc[storIdx,'q_mvar'].item()\n",
    "\n",
    "        S_stor = mt.sqrt(P_stor**2+Q_stor**2)\n",
    "\n",
    "        SOC_ini_var = (S_stor*period_duration.item())/Storlist[i].max_ener\n",
    "\n",
    "        SOC_f = SOC_ini_var + Storlist[i].ini_SOC\n",
    "\n",
    "        if SOC_f < Storlist[i].min_SOC: # SOC lower limit verification \n",
    "\n",
    "            SOC_max_var = Storlist[i].ini_SOC - Storlist[i].min_SOC\n",
    "            net.storage.loc[storIdx,'p_mw'] = (SOC_max_var*Storlist[i].max_ener/period_duration.item())*Storlist[i].pn/1000\n",
    "\n",
    "        elif SOC_f > Storlist[i].max_SOC: # SOC upper limit verification\n",
    "\n",
    "            SOC_max_var = Storlist[i].max_SOC - Storlist[i].ini_SOC\n",
    "            net.storage.loc[storIdx,'p_mw'] = (SOC_max_var*Storlist[i].max_ener/period_duration.item())*Storlist[i].pn/1000\n",
    "\n",
    "        Storlist[i].ini_SOC = SOC_f\n",
    "\n",
    "        net.storage.loc[storIdx,'soc_percent'] = Storlist[i].ini_SOC        \n",
    "\n",
    "        aux_Matrixsoc[i,t] = net.storage.loc[storIdx,'soc_percent'].item()\n",
    "\n",
    "    return net,aux_Matrixsoc\n",
    "\n",
    "def convEst(net,Nodelist,Convlist,bus_voltages,bus_voltages_angles,base_power,tol_max_perc): # In progress\n",
    "\n",
    "    for i in range(0,len(Convlist)):\n",
    "\n",
    "        if Convlist[i].convtype == 'Converter DC/DC':\n",
    "\n",
    "            convIdx = net.trafo.index[('Converter DC/DC ' + str(Convlist[i].conv_id))==net.trafo['name']].item()\n",
    "\n",
    "            if (bus_voltages[buses[Convlist[i].from_bus]]-bus_voltages[buses[Convlist[i].to_bus]]) == 0:\n",
    "\n",
    "                Req_corr = 0.005*(Convlist[i].pn/base_power.item())\n",
    "\n",
    "            else:\n",
    "\n",
    "                Pconv = (bus_voltages[buses[Convlist[i].from_bus]]-bus_voltages[buses[Convlist[i].to_bus]])**2/net.trafo.loc[convIdx,'vk_percent']\n",
    "                \n",
    "                aux_PF = Convlist[i].eff_table['Efficiency Curve'].dropna().unique()\n",
    "                auxidxEC = Convlist[i].eff_table.columns.get_loc('Efficiency Curve')\n",
    "                aux_PO = Convlist[i].eff_table.columns[auxidxEC+1:].values\n",
    "                aux_eff = Convlist[i].eff_table.iloc[:,auxidxEC+1:]\n",
    "\n",
    "                points = [(fp, p) for fp in aux_PF for p in aux_PO]\n",
    "                values = aux_eff.values.flatten()\n",
    "\n",
    "                fp_novo, p_novo = 1, Pconv\n",
    "\n",
    "                eff = griddata(points, values, (fp_novo, p_novo), method='linear')\n",
    "\n",
    "                Req = bus_voltages[buses[Convlist[i].from_bus]]*(bus_voltages[buses[Convlist[i].from_bus]]-bus_voltages[buses[Convlist[i].to_bus]])/(net.res_bus.p_mw[buses[Convlist[i].from_bus]].item()/(base_power.item()*1000))\n",
    "\n",
    "                if abs(Req-(((bus_voltages[buses[Convlist[i].from_bus]].item()-bus_voltages[buses[Convlist[i].to_bus]].item())**2)/((1-eff)*(net.res_bus.p_mw[buses[Convlist[i].from_bus]]/(base_power.item()*1000)))))>tol_max_perc:\n",
    "\n",
    "                    Req_corr = (((bus_voltages[buses[Convlist[i].from_bus]].item()-bus_voltages[buses[Convlist[i].to_bus]].item())**2)/((1-eff)*(net.res_bus.p_mw[buses[Convlist[i].from_bus]]/(base_power.item()*1000)))) # In grid pu\n",
    "                \n",
    "                else:\n",
    "\n",
    "                    Req_corr = Req\n",
    "\n",
    "            Req_corr = Req_corr * (base_power.item()/Convlist[i].pn)\n",
    "\n",
    "            net.trafo.loc[convIdx,'vk_percent'] = 1.001*Req_corr\n",
    "\n",
    "            net.trafo.loc[convIdx,'vkr_percent'] = Req_corr\n",
    "        \n",
    "        if Convlist[i].convtype == 'Converter AC/DC': \n",
    "\n",
    "            conv_trafoIdx = net.trafo.index[('Converter AC/DC ' + str(Convlist[i].conv_id))==net.trafo['name']].item() # Index of the converter in transformer list in net data\n",
    "\n",
    "            conv_wardIdx = net.ward.index[('Converter AC/DC ' + str(Convlist[i].conv_id))==net.ward['name']].item() # Index of the converter in ward list in net data\n",
    "            \n",
    "            bus_comp_volt = []\n",
    "\n",
    "            k = 1.01    # Constant just to make the ratio between the resistance and the reactance of the transformer model in order to allow the convergence of the power flow\n",
    "\n",
    "            for j in range(0,len(bus_voltages)):\n",
    "                bus_comp_volt.append(complex(bus_voltages[j]*mt.cos(bus_voltages_angles[j]),bus_voltages[j]*mt.sin(bus_voltages_angles[j])))    # Obtaining the complex voltages of each bus\n",
    "\n",
    "            if all(i == 1 for i in bus_comp_volt):\n",
    "                \n",
    "                Ractual = 1e-3                                  # Actual resistance (pu grid)\n",
    "                absZactual = k*Ractual\n",
    "                Xactual = mt.sqrt(absZactual**2-Ractual**2)     # Actual reactance (pu grid)\n",
    "                Req_corr = Ractual\n",
    "                Xeq_corr = Xactual\n",
    "                ps_conv = 0\n",
    "\n",
    "            else:\n",
    "                \n",
    "                # Actualized iteration converter parameters\n",
    "                 \n",
    "                absZactual = (net.trafo.loc[conv_trafoIdx,'vk_percent']/100) * (base_power.item()/(Convlist[i].pn)) / (1e3)     # Actual impedance (pu grid)    (1e3 due to Pandapower works in MVA instead of kVA)\n",
    "                Ractual = (net.trafo.loc[conv_trafoIdx,'vkr_percent']/100) * (base_power.item()/(Convlist[i].pn)) / (1e3)       # Actual resistance (pu grid)   \n",
    "                Xactual = mt.sqrt(absZactual**2-Ractual**2)                                                                     # Actual reactance (pu grid)\n",
    "                Zactual = complex(Ractual,Xactual)\n",
    "                Yactual = 1/Zactual\n",
    "                ps_actual = net.trafo.loc[conv_trafoIdx,'shift_degree'] * (mt.pi/180)                                           # Phase-shift parameter in transformer converter equivalent model (Already save in a variable in case that if needed)\n",
    "                m = 1                                                                                                           # We assume for now a constant ratio between voltages \n",
    "                m_complex_actual = complex(m*np.cos(ps_actual),m*np.sin(ps_actual))\n",
    "                Qs_ward_actual = net.ward.loc[conv_wardIdx,'qs_mvar']\n",
    "\n",
    "                # Verification of lv, hv, AC and DC buses and their voltages of the 2 elements converter model\n",
    "\n",
    "                def verf_buses(net,conv_trafoIdx,conv_wardIdx):\n",
    "\n",
    "                    hv_bus = net.trafo.loc[conv_trafoIdx,'hv_bus'].item()\n",
    "                    lv_bus = net.trafo.loc[conv_trafoIdx,'lv_bus'].item()\n",
    "                    AC_bus = net.ward.loc[conv_wardIdx,'bus'].item()\n",
    "                    \n",
    "                    if AC_bus==hv_bus:\n",
    "                        DC_bus = lv_bus\n",
    "                    else:\n",
    "                        DC_bus = hv_bus\n",
    "\n",
    "                    return hv_bus, lv_bus, AC_bus, DC_bus\n",
    "\n",
    "                # Function for computation of Power Flows of Pandapower Transformer model (HV side to LV side and vice-versa)\n",
    "\n",
    "                def basic_power_flow_eqs(hv_bus,lv_bus,V_hv,V_lv,Yactual,m_complex_actual):\n",
    "\n",
    "                    S_hv_lv = (Yactual.conjugate()/(abs(m_complex_actual)**2))*(abs(V_hv)**2-m_complex_actual.conjugate()*V_hv*V_lv.conjugate())                    # Complex Power flowing from HV side to LV side \n",
    "                    S_lv_hv = (Yactual.conjugate()/(abs(m_complex_actual)**2))*((abs(V_lv)**2)*(abs(m_complex_actual)**2)-m_complex_actual*V_lv*V_hv.conjugate())   # Complex Power flowing from LV side to HV side\n",
    "\n",
    "                    S_losses = S_hv_lv + S_lv_hv\n",
    "                    P_losses = S_losses.real\n",
    "                    Q_losses = S_losses.imag\n",
    "\n",
    "                    return S_hv_lv, S_lv_hv, P_losses, Q_losses\n",
    "\n",
    "                # Function for verification of power flow (if power flows from AC to DC or vice-versa)\n",
    "\n",
    "                def verf_power_flow(S_hv_lv,S_lv_hv):\n",
    "\n",
    "                    if abs(S_hv_lv.real) > abs(S_lv_hv.real):\n",
    "                        Pout = S_lv_hv.real\n",
    "                        Qout = S_lv_hv.imag\n",
    "                        Sout = complex(Pout,Qout)\n",
    "                        FPout = mt.cos(abs(Pout/abs(Sout)))\n",
    "                        Pin = S_hv_lv.real\n",
    "                        Qin = S_hv_lv.imag\n",
    "                        Sin = complex(Pin,Qin)\n",
    "                        FPin = mt.cos(abs(Pin/abs(Sin)))\n",
    "                        eff_actual = abs(Pout/Pin)\n",
    "\n",
    "                    else:\n",
    "                        Pout = S_hv_lv.real\n",
    "                        Qout = S_hv_lv.imag\n",
    "                        Sout = complex(Pout,Qout)\n",
    "                        FPout = mt.cos(abs(Pout/abs(Sout)))\n",
    "                        Pin = S_lv_hv.real\n",
    "                        Qin = S_lv_hv.imag\n",
    "                        Sin = complex(Pin,Qin)\n",
    "                        FPin = mt.cos(abs(Pin/abs(Sin)))\n",
    "                        eff_actual = abs(Pout/Pin)\n",
    "                        \n",
    "                    return Sout,FPout,Sin\n",
    "\n",
    "                # Estimation of efficiency from converter data\n",
    "\n",
    "                def effest(Pout,FPOut,aux_PF,aux_PO,aux_eff):\n",
    "\n",
    "                    points = [(fp, p) for fp in aux_PF for p in aux_PO]\n",
    "                    values = aux_eff.values.flatten()\n",
    "\n",
    "                    fp_novo, p_novo = FPOut, Pout\n",
    "\n",
    "                    eff = griddata(points, values, (fp_novo, p_novo), method='linear')\n",
    "\n",
    "                    return eff\n",
    "\n",
    "                ######## DECISION/CONTROL MODULE OF CONVERTER ########\n",
    "\n",
    "                # Preliminary computations \n",
    "\n",
    "                hv_bus, lv_bus, AC_bus, DC_bus = verf_buses(net,conv_trafoIdx,conv_wardIdx)\n",
    "                V_hv = bus_comp_volt[hv_bus]\n",
    "                V_lv = bus_comp_volt[lv_bus]\n",
    "                V_AC = bus_comp_volt[AC_bus]\n",
    "                V_DC = bus_comp_volt[DC_bus]\n",
    "                S_hv_lv, S_lv_hv, P_losses, Q_losses = basic_power_flow_eqs(hv_bus,lv_bus,V_hv,V_lv,Yactual,m_complex_actual)\n",
    "                Sout, FPout, Sin = verf_power_flow(S_hv_lv,S_lv_hv)\n",
    "\n",
    "                # Auxiliary vectors and indexes for the determination of the point of efficiency of the converter\n",
    "\n",
    "                aux_PF = Convlist[i].eff_table['Efficiency Curve'].dropna().unique()\n",
    "                auxidxEC = Convlist[i].eff_table.columns.get_loc('Efficiency Curve')\n",
    "                aux_PO = Convlist[i].eff_table.columns[auxidxEC+1:].values\n",
    "                aux_eff = Convlist[i].eff_table.iloc[:,auxidxEC+1:]\n",
    "\n",
    "                # Point of operation interpolated from efficiency table from converter data\n",
    "\n",
    "                Pout = abs(Sout.real)\n",
    "                Pin = abs(Sin.real)\n",
    "                eff_actual = effest(Pout*(base_power.item()/(Convlist[i].pn)),FPout,aux_PF,aux_PO,aux_eff)\n",
    "\n",
    "                # Control module\n",
    "\n",
    "                # Verification functions (In progress)\n",
    "\n",
    "                if eff_actual != Pout/Pin:  # Function to verify the efficiency of pandapower vs efficiency table \n",
    "                    \n",
    "                    Gt = (Pin * (1-eff_actual))/(((abs(V_AC)**2+(abs(m_complex_actual)**2)*abs(V_DC)**2)/abs(m_complex_actual)**2)+(((2*abs(V_AC)*abs(V_DC))/abs(m_complex_actual))*np.cos(np.angle(V_AC)-np.angle(V_DC)-np.angle(m_complex_actual))))\n",
    "\n",
    "                    print(Pin * (1-eff_actual),P_losses,((abs(V_AC)**2+(abs(m_complex_actual)**2)*abs(V_DC)**2)/abs(m_complex_actual)**2)+(((2*abs(V_AC)*abs(V_DC))/abs(m_complex_actual))*np.cos(np.angle(V_AC)-np.angle(V_DC)-np.angle(m_complex_actual))))\n",
    "\n",
    "                    Ractual = 1/(Gt*(k**2))\n",
    "\n",
    "                    Xactual = Ractual*mt.sqrt(k**2-1)\n",
    "\n",
    "                if abs(S_hv_lv.imag) > abs(S_lv_hv.imag) and abs(S_hv_lv.imag)>1e-6:    # Function to verify the reactive power flow in branch of transformer model and compensate it   \n",
    "\n",
    "                    Qs_ward_actual = S_hv_lv.imag\n",
    "\n",
    "                elif abs(S_lv_hv.imag) > abs(S_hv_lv.imag) and abs(S_lv_hv.imag)>1e-6:\n",
    "\n",
    "                    Qs_ward_actual = S_lv_hv.imag\n",
    "\n",
    "                else:\n",
    "\n",
    "                    Qs_ward_actual = Qs_ward_actual\n",
    "                     \n",
    "                # Control function (To implement a function that allows the output voltage or power to be controlled by a reference value for exemple a imposed voltage output or imposed power output) (In progress)\n",
    "                   \n",
    "\n",
    "                Req_corr = Ractual\n",
    "                Xeq_corr = Xactual\n",
    "                ps_conv = ps_actual * (180/mt.pi)\n",
    "                Qs_ward_corr = Qs_ward_actual\n",
    "\n",
    "                print('\\nR = ',Req_corr,'\\nX = ',Xeq_corr,'\\nEficiencia atual= ',eff_actual)\n",
    "\n",
    "\n",
    "            Req_corr = Req_corr * (Convlist[i].pn/base_power.item())\n",
    "\n",
    "            Xeq_corr = Xeq_corr * (Convlist[i].pn/base_power.item())\n",
    "\n",
    "            net.trafo.loc[conv_trafoIdx,'vk_percent'] = mt.sqrt(Req_corr**2+Xeq_corr**2) * 100\n",
    "\n",
    "            net.trafo.loc[conv_trafoIdx,'vkr_percent'] = Req_corr * 100\n",
    "\n",
    "            net.trafo.loc[conv_trafoIdx,'shift_degree'] = ps_conv\n",
    "\n",
    "    return net\n",
    "\n",
    "\n",
    "def runpptd_cusm(net,buses,Nodelist,Linelist,Trafolist,Convlist,Loadlist,Genlist,Storlist,base_power):\n",
    "\n",
    "    time_list_chosen = timeQuest(time_list) # List of time instant(s) chosen by the user \n",
    "\n",
    "    aux_forloop = False # Auxiliary variable to know if the for loop of T is running as the first time (assumes a PF start as a 'flat start') or have previous time values\n",
    "                        # Since the previous solution has already set the power values to loads, in the run of the next time iteration is not needed to prepare the values using the new VP droop curves\n",
    "\n",
    "    tol_max_perc = 1e-8 # Maximum tolerance in percentage\n",
    "\n",
    "    n_iter_max = 10000 # Maximum number of iterations assumed\n",
    "\n",
    "    n_iter_PF = 10 # Maximum number of iteration of Pandapower PF \n",
    "\n",
    "    erro = 1 # Error\n",
    "\n",
    "    bus_voltages_old_print = [np.zeros([len(Nodelist),n_iter_max]) for i in time_list_chosen] # Just for printing bus voltages behaviors with the evolution of the iterations\n",
    "\n",
    "    load_p_old_print = [np.zeros([len(Loadlist),n_iter_max]) for i in time_list_chosen]         # Just for printing how the power in each load behaves in each iteration and in each time iteration\n",
    "\n",
    "    SOC_old_print = [np.zeros([len(Storlist),n_iter_max]) for i in time_list_chosen]         # Just for printing how the SOC in each BESS behaves in each iteration and in each time iteration\n",
    "\n",
    "    aux_Matrix = np.ones([len(Nodelist),n_iter_max])    # Auxiliary matrix for bus voltages \n",
    "\n",
    "    aux_Matrixp = np.ones([len(Loadlist),n_iter_max])   # Auxiliary matrix for loads power\n",
    "\n",
    "    aux_Matrixsoc = np.ones([len(Storlist),n_iter_max]) # Auxiliary matrix for SOC \n",
    "\n",
    "    start_time_TD = tm.time()\n",
    "\n",
    "    for T in time_list_chosen:\n",
    "\n",
    "        T = int(T)\n",
    "\n",
    "        # Run power flow analysis\n",
    "        print(\">Run power flow analysis\")\n",
    "\n",
    "        repetir = True # Condition to continue while loop until the PF analysis is not completed\n",
    "\n",
    "        t = 0 # Iteration counter\n",
    "\n",
    "        # Droop preparation\n",
    "\n",
    "        bus_voltages_old = np.zeros(len(Nodelist))         # Auxiliary array to memorize the results for comparison them between iterations \n",
    "        bus_power_injections_old = np.zeros(len(Nodelist))         # Auxiliary array to memorize the results for comparison them between iterations \n",
    "\n",
    "        for i in range(0,len(Nodelist)):\n",
    "        \n",
    "            aux_Matrix[buses[Nodelist[i].node_number.item()],t] = 1               # Define the values of the first iteration for printing \n",
    "        \n",
    "        if aux_forloop is False:                        # Droop adjustment or converters state estimation for initial time step \n",
    "\n",
    "            # Converter equipment\n",
    "\n",
    "            net = convEst(net,Nodelist,Convlist,np.ones(len(Nodelist)),np.zeros(len(Nodelist)),base_power,tol_max_perc)\n",
    "            \n",
    "            # Load equipment \n",
    "\n",
    "            [net,aux_Matrixp] = loaddroopReadjustment(net,Loadlist,np.ones(len(Nodelist)),base_power,aux_Matrixp,T,t)\n",
    "\n",
    "            # Generator equipment\n",
    "\n",
    "            net = gendroopReadjustment(net,Genlist,np.ones(len(Nodelist)),base_power,T)\n",
    "\n",
    "            # Storage equipment\n",
    "\n",
    "            [net,aux_Matrixsoc] = storControl(net,Storlist,np.ones(len(Nodelist)),base_power,period_duration,aux_Matrixsoc,T,t)\n",
    "\n",
    "            aux_forloop = True\n",
    "\n",
    "        start_time = tm.time()\n",
    "\n",
    "        while (repetir is True):\n",
    "\n",
    "            ###### PowerFlow analysis procedure ######\n",
    "\n",
    "            try:\n",
    "                pp.runpp(net,max_iteration=n_iter_PF)\n",
    "                print(net.res_bus)\n",
    "            except pp.LoadflowNotConverged:\n",
    "                print(\"Power Flow did not converge\")\n",
    "\n",
    "            bus_indices = net.res_bus.index + 1     # adjust to start from 1\n",
    "            bus_voltages = net.res_bus.vm_pu\n",
    "            bus_voltages_angles = net.res_bus.va_degree * (mt.pi/180)\n",
    "            bus_power_injections = net.res_bus.p_mw\n",
    "            \n",
    "            desvio = abs(bus_voltages-bus_voltages_old)/bus_voltages + abs(bus_power_injections-bus_power_injections_old)/bus_power_injections # Deviation measurement of voltage profile of each iteration and active power injections (in percentage) \n",
    "\n",
    "            erro = (max(desvio))*100 # Error computation of the higher deviation in percentage\n",
    "\n",
    "            if (erro < tol_max_perc) and t!=0:               # Comparison between iterations resorting to error\n",
    "                end_time = tm.time()\n",
    "                elapsed_time = end_time - start_time\n",
    "                print(\"Power Flow solution obtained with an error of\",erro,\"% at\",t+1,\"iterations in\",elapsed_time*1000,\"ms for time instant\",T)\n",
    "                repetir = False\n",
    "\n",
    "            elif (t>=n_iter_max):\n",
    "                print(\"Power Flow solution not obtained due to maximum number of iterations being reached with an error of\",erro,\"%\")\n",
    "                repetir = False\n",
    "\n",
    "            bus_voltages_old = bus_voltages.values\n",
    "            bus_power_injections_old = bus_power_injections.values\n",
    "\n",
    "            ###### Control Functions - Converter Equivalent Resistance Estimation and Droop Readjusment ######\n",
    "            \n",
    "            # Converter equipment\n",
    "\n",
    "            net = convEst(net,Nodelist,Convlist,bus_voltages,bus_voltages_angles,base_power,tol_max_perc)\n",
    "\n",
    "            # Load equipment\n",
    "\n",
    "            [net,aux_Matrixp] = loaddroopReadjustment(net,Loadlist,bus_voltages,base_power,aux_Matrixp,T,t)\n",
    "                \n",
    "            # Generator equipment\n",
    "\n",
    "            net = gendroopReadjustment(net,Genlist,bus_voltages,base_power,T)\n",
    "\n",
    "            # Storage equipment\n",
    "\n",
    "            [net,aux_Matrixsoc] = storControl(net,Storlist,bus_voltages,base_power,period_duration,aux_Matrixsoc,T,t)\n",
    "\n",
    "            t=t+1 # Iteration counter\n",
    "\n",
    "            # Define the values for the atual iteration for printing\n",
    "\n",
    "            for i in range(0,len(Nodelist)):\n",
    "\n",
    "                aux_Matrix[buses[Nodelist[i].node_number.item()],t] = bus_voltages_old[buses[Nodelist[i].node_number.item()]]        # Define the values of the first iteration for printing \n",
    "\n",
    "            # Remove the remaining elements that does not count for printing purposes\n",
    "\n",
    "            if (repetir == False):\n",
    "                \n",
    "                aux_Matrix = np.delete(aux_Matrix,range(t,n_iter_max),1) \n",
    "                aux_Matrixp = np.delete(aux_Matrixp,range(t,n_iter_max),1)\n",
    "                aux_Matrixsoc = np.delete(aux_Matrixsoc,range(t,n_iter_max),1)\n",
    "        \n",
    "        if len(time_list_chosen)>1:\n",
    "\n",
    "            bus_voltages_old_print[T-1] = aux_Matrix      # Since the index might be not exact with the time instant for exemple (time instant chosen of 2 and 3, so first index is not 2), this must be improved\n",
    "            load_p_old_print[T-1] = aux_Matrixp\n",
    "            SOC_old_print[T-1] = aux_Matrixsoc\n",
    "\n",
    "        else:\n",
    "\n",
    "            bus_voltages_old_print = aux_Matrix\n",
    "            load_p_old_print = aux_Matrixp\n",
    "            SOC_old_print = aux_Matrixsoc\n",
    "\n",
    "        del aux_Matrix,aux_Matrixp,aux_Matrixsoc\n",
    "\n",
    "        aux_Matrix = np.zeros([len(Nodelist),n_iter_max])       # It is needed to resize the aux_Matrix since the loop is going to go to another time instant\n",
    "        aux_Matrixp = np.zeros([len(Loadlist),n_iter_max])      # It is needed to resize the aux_Matrixp since the loop is going to go to another time instant\n",
    "        aux_Matrixsoc = np.zeros([len(Storlist),n_iter_max])    # It is needed to resize the aux_Matrixsoc since the loop is going to go to another time instant\n",
    "\n",
    "    end_time_TD = tm.time()\n",
    "    elapsed_time_TD = end_time_TD - start_time_TD\n",
    "\n",
    "    print(\"Power Flow solution obtained in\",elapsed_time_TD*1000,\"ms for time interval between\",min(time_list_chosen),'and',max(time_list_chosen))\n",
    "\n",
    "    return net,bus_voltages_old_print,load_p_old_print,SOC_old_print,bus_indices,bus_voltages,bus_voltages_angles\n",
    "\n",
    "[net,bus_voltages_old_print,load_p_old_print,SOC_old_print,bus_indices,bus_voltages,bus_voltages_angles] = runpptd_cusm(net,buses,Nodelist,Linelist,Trafolist,Convlist,Loadlist,Genlist,Storlist,base_power)\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Print results\n",
    "print(\">Print results\")\n",
    "print(net.res_bus)\n",
    "print(net.res_line)\n",
    "print(net.res_storage)\n",
    "print(net.res_trafo)\n",
    "# print(\"Bus voltages (pu):\", net.res_bus.vm_pu)\n",
    "# print(\"Line loading (A):\", net.res_line.loading_percent)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot voltage profiles\n",
    "\n",
    "def printVoltProf(bus_indices,bus_voltages_old_print):\n",
    "    print(\">Plot voltage profiles\")\n",
    "\n",
    "    fig, axs = plt.subplots(len(bus_voltages_old_print),1,figsize=(14, 30),layout='constrained')\n",
    "    title_names = []\n",
    "\n",
    "    for i in range(0,len(bus_voltages_old_print)):       \n",
    "        aux_data = bus_voltages_old_print[i]\n",
    "        if aux_data.shape[1]>1:\n",
    "            axs[i].plot(bus_indices,aux_data[:,-1])\n",
    "        else:\n",
    "            axs[i].plot(bus_indices,aux_data[:])\n",
    "        title_names = ['Time Instant = ' + str(i + 1)]\n",
    "        axs[i].set_title(title_names)\n",
    "        axs[i].set_ylabel(\"Voltage (pu)\")\n",
    "        axs[i].grid(True)\n",
    "\n",
    "\n",
    "    #plt.figure(figsize=(8.5, 2.25))\n",
    "    #plt.plot(bus_indices, bus_voltages_old_print, marker='o', linestyle='-', color='b')\n",
    "    #plt.title(\"Voltage Profiles\")\n",
    "    #plt.xlabel(\"Bus Index\")\n",
    "    #plt.ylabel(\"Voltage (pu)\")\n",
    "    #plt.grid(True)\n",
    "    # plt.ylim(0.9, 1.1)\n",
    "    # plt.ylim(0.98, 1.06)\n",
    "    #plt.xticks(bus_indices)\n",
    "    plt.show()\n",
    "    return\n",
    "\n",
    "printVoltProf(bus_indices,bus_voltages_old_print)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot bus voltage evolution in iterations by each bus \n",
    "\n",
    "def printVoltIter(bus_voltages_old_print):\n",
    "    print(\">Plot bus voltages vs iterations evolution\")\n",
    "\n",
    "    if type(bus_voltages_old_print) is list:\n",
    "        \n",
    "        auxPrint = 0\n",
    "        auxIdx = []\n",
    "        \n",
    "        for i in range(0,len(bus_voltages_old_print)):\n",
    "            \n",
    "            if bus_voltages_old_print[i].shape[1] > 1:\n",
    "                \n",
    "                auxPrint += 1\n",
    "                auxIdx.append(i)\n",
    "\n",
    "\n",
    "        fig, axs = plt.subplots(auxPrint,1,figsize=(14, 30),layout='constrained')\n",
    "        title_names = []\n",
    "\n",
    "        for i in range(0,auxPrint):\n",
    "                \n",
    "                axs[i].plot(bus_voltages_old_print[auxIdx[i]].transpose())\n",
    "                title_names = ['Time Instant = ' + str(auxIdx[i] + 1)]\n",
    "                axs[i].set_title(title_names)\n",
    "                axs[i].set_ylabel(\"Voltage (pu)\")\n",
    "                axs[i].set_xlim(0, bus_voltages_old_print[auxIdx[i]].shape[1]-1)\n",
    "                axs[i].set_xticks(range(1,bus_voltages_old_print[auxIdx[i]].shape[1]))\n",
    "                axs[i].grid(True)\n",
    "        \n",
    "        plt.xlabel(\"Iterations\")\n",
    "        plt.show\n",
    "\n",
    "    else:\n",
    "\n",
    "        plt.plot(bus_voltages_old_print.transpose())\n",
    "        plt.show\n",
    "\n",
    "    return\n",
    "\n",
    "def printPowerIter(load_p_old_print):\n",
    "    \n",
    "    print(\">Plot power components vs iterations evolution\")\n",
    "\n",
    "    if type(load_p_old_print) is list:\n",
    "        \n",
    "        auxPrint = 0\n",
    "        auxIdx = []\n",
    "        \n",
    "        for i in range(0,len(load_p_old_print)):\n",
    "            \n",
    "            if load_p_old_print[i].shape[1] > 1:\n",
    "                \n",
    "                auxPrint += 1\n",
    "                auxIdx.append(i)\n",
    "\n",
    "\n",
    "        fig, axs = plt.subplots(auxPrint,1,figsize=(14, 30),layout='constrained')\n",
    "        title_names = []\n",
    "\n",
    "        for i in range(0,auxPrint):\n",
    "                \n",
    "                axs[i].plot(load_p_old_print[auxIdx[i]].transpose())\n",
    "                title_names = ['Time Instant = ' + str(auxIdx[i] + 1)]\n",
    "                axs[i].set_title(title_names)\n",
    "                axs[i].set_ylabel(\"Power (pu)\")\n",
    "                axs[i].set_xlim(0, load_p_old_print[auxIdx[i]].shape[1]-1)\n",
    "                axs[i].set_xticks(range(1,load_p_old_print[auxIdx[i]].shape[1]))\n",
    "                axs[i].grid(True)\n",
    "        \n",
    "        plt.xlabel(\"Iterations\")\n",
    "        plt.show\n",
    "\n",
    "    else:\n",
    "\n",
    "        plt.plot(load_p_old_print.transpose())\n",
    "        plt.show\n",
    "\n",
    "    return\n",
    "\n",
    "def printSOCIter(SOC_old_print):\n",
    "    \n",
    "    print(\">Plot SOC vs iterations evolution\")\n",
    "\n",
    "    if type(SOC_old_print) is list:\n",
    "        \n",
    "        auxPrint = 0\n",
    "        auxIdx = []\n",
    "        \n",
    "        for i in range(0,len(SOC_old_print)):\n",
    "            \n",
    "            if SOC_old_print[i].shape[1] > 1:\n",
    "                \n",
    "                auxPrint += 1\n",
    "                auxIdx.append(i)\n",
    "\n",
    "\n",
    "        fig, axs = plt.subplots(auxPrint,1,figsize=(14, 30),layout='constrained')\n",
    "        title_names = []\n",
    "\n",
    "        for i in range(0,auxPrint):\n",
    "                \n",
    "                axs[i].plot(SOC_old_print[auxIdx[i]].transpose())\n",
    "                title_names = ['Time Instant = ' + str(auxIdx[i] + 1)]\n",
    "                axs[i].set_title(title_names)\n",
    "                axs[i].set_ylabel(\"SOC (%)\")\n",
    "                axs[i].set_xlim(0, SOC_old_print[auxIdx[i]].shape[1]-1)\n",
    "                axs[i].set_xticks(range(1,SOC_old_print[auxIdx[i]].shape[1]))\n",
    "                axs[i].grid(True)\n",
    "        \n",
    "        plt.xlabel(\"Iterations\")\n",
    "        plt.show\n",
    "\n",
    "    else:\n",
    "\n",
    "        plt.plot(SOC_old_print.transpose())\n",
    "        plt.show\n",
    "\n",
    "    return\n",
    "\n",
    "printVoltIter(bus_voltages_old_print)\n",
    "\n",
    "printPowerIter(load_p_old_print)\n",
    "\n",
    "printSOCIter(SOC_old_print)\n",
    "\n",
    "''' Code written that could help\n",
    "\n",
    "fig, axs = plt.subplots(len(drp_idx_aux),1,figsize=(14, 30),layout='constrained')\n",
    "#fig.title(\"Droop Evolution\")\n",
    "labels = []\n",
    "for i,_ in enumerate(droop_data.Node_j):\n",
    "    axs[i].plot(range(0,t), bus_voltages_old_print[:,i], marker='o', linestyle='-')\n",
    "    labels = [net.bus['name'].iloc[int(drp_idx_aux[i])]+' - PF']\n",
    "    axs[i].set_ylabel(\"Voltage (pu)\")\n",
    "    axs[i].legend(labels)\n",
    "    axs[i].grid(True)\n",
    "    axs[i].set_xlim(0, t)\n",
    "    #axs[i].set_ylim(bus_voltages_old_print[t-1,i]*0.995, bus_voltages_old_print[t-1,i]*1.005)\n",
    "    #plt.plot(range(1,t+1), volt_aux.values[i,:], marker='+', linestyle=':')    \n",
    "    #labels.append([net.bus['name'].iloc[int(drp_idx_aux[i])]+' - Power Voltage Droop Curve'])\n",
    "plt.xlabel(\"Iterations\")\n",
    "\n",
    "labels = []\n",
    "#plt.ylim(0.925, 1.075)\n",
    "plt.show()\n",
    "\n",
    "print(\">Plot droop evolution\")\n",
    "\n",
    "plt.figure(figsize=(14, 6))\n",
    "plt.title(\"Droop Evolution\")\n",
    "labels = []\n",
    "for i,_ in enumerate(droop_data.Node_j):\n",
    "    plt.plot(load_p_old_print[:,i],bus_voltages_old_print[:,i], marker='o', linestyle='-')\n",
    "    labels.append([net.bus['name'].iloc[int(drp_idx_aux[i])]+' - PF'])\n",
    "    plt.plot(p_aux.values[i,:], volt_aux.values[i,:], marker='+', linestyle=':')    \n",
    "    labels.append([net.bus['name'].iloc[int(drp_idx_aux[i])]+' - Power Voltage Droop Curve'])\n",
    "plt.legend(labels)\n",
    "plt.xlabel(\"Power (pu)\")\n",
    "plt.ylabel(\"Voltage (pu)\")\n",
    "plt.grid(True)\n",
    "#plt.ylim(0.925, 1.075)\n",
    "#plt.xlim(-0.05, 2.1)\n",
    "# plt.xticks(bus_indices)\n",
    "plt.show()\n",
    "\n",
    "'''"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Computation of Environmental KPIs\n",
    "\n",
    "# PUE - Power Usage Effectiveness (by load component)\n",
    "\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "base",
   "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.12.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}