{ "cells": [ { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " ----------------------------------------------------------------\n", " APMonitor, Version 1.0.3\n", " APMonitor Optimization Suite\n", " ----------------------------------------------------------------\n", " \n", " \n", " --------- APM Model Size ------------\n", " Each time step contains\n", " Objects : 97\n", " Constants : 0\n", " Variables : 312\n", " Intermediates: 101\n", " Connections : 291\n", " Equations : 306\n", " Residuals : 205\n", " \n", " Number of state variables: 500\n", " Number of total equations: - 398\n", " Number of slack variables: - 7\n", " ---------------------------------------\n", " Degrees of freedom : 95\n", " \n", " ----------------------------------------------\n", " Model Parameter Estimation with APOPT Solver\n", " ----------------------------------------------\n", "Iter: 1 I: 0 Tm: 1.47 NLPi: 55 Dpth: 0 Lvs: 2 Obj: 0.00E+00 Gap: NaN\n", "Iter: 2 I: -9 Tm: 273.42 NLPi:10001 Dpth: 1 Lvs: 1 Obj: 0.00E+00 Gap: NaN\n", "--Integer Solution: 4.64E+03 Lowest Leaf: 4.64E+03 Gap: 0.00E+00\n", "Iter: 3 I: 0 Tm: 0.10 NLPi: 4 Dpth: 1 Lvs: 1 Obj: 4.64E+03 Gap: 0.00E+00\n", " Successful solution\n", " \n", " ---------------------------------------------------\n", " Solver : APOPT (v1.0)\n", " Solution time : 275.001799999998 sec\n", " Objective : 4642.69974872768 \n", " Successful solution\n", " ---------------------------------------------------\n", " \n", "\n", "\n", "Optimized design:\n", "\tFreq [MHz]: 13.560\n", "\tDo_1 [mm]: 200.000\n", "\tDo_2 [mm]: 200.000\n", "\tGap_Distance [mm]: 2.732\n", "\tPitch_1 [mm]: 8.730\n", "\tPitch_2 [mm]: 4.752\n", "\tTurns: 4.000\n", "\tWidth_1 [mm]: 8.730\n", "\tWidth_2 [mm]: 11.693\n", "\tDi_1 [mm]: 80.564\n", "\tDi_2 [mm]: 80.564\n", "\tThickness_1 [mm]: 0.500\n", "\tThickness_2 [mm]: 0.500\n", "Maximized Q_factor: -4642.6997487\n", "\u001b[1m1/1\u001b[0m \u001b[32m━━━━━━━━━━━━━━━━━━━━\u001b[0m\u001b[37m\u001b[0m \u001b[1m0s\u001b[0m 31ms/step\n", "Maximized Q_factor: 1017.166\n", "Gekko Solvetime: 275.0018 sec\n", "Gekko Solvetime: 00:04:35\n", "L_coil= 2.483 uH\n", "C_coil= 53.258 pF\n", "Resonant Freq = 13.841 MHz\n" ] } ], "source": [ "import numpy as np\n", "import pandas as pd\n", "from gekko import GEKKO\n", "from gekko.ML import Gekko_NN_TF, CustomMinMaxGekkoScaler\n", "from tensorflow.keras.models import load_model\n", "import warnings\n", "warnings.filterwarnings(\"ignore\")\n", "\n", "# Q-factor prediction from trained TF model\n", "model_path = '/Users/UW/Summer 2024/Thesis/04_ML/Model/v3/best_model_90_10_withInit.keras'\n", "model = load_model(model_path)\n", "\n", "# CSV file reading\n", "data_path = '/Users//UW/Summer 2024/Thesis/04_ML/Total_data/Total_data.csv'\n", "data = pd.read_csv(data_path)\n", "\n", "# Specify feature columns and target column\n", "target = ['Q_factor']\n", "features = [col for col in data.columns if col not in target]\n", "\n", "# Find max/min of scaler for Gekko_NN_TF parameters\n", "s = CustomMinMaxGekkoScaler(data, features, target)\n", "mma = s.minMaxValues()\n", "\n", "# Initialize GEKKO\n", "m = GEKKO(remote=False)\n", "\n", "# Design Constants\n", "# epsilon = m.Const(value=8.854e-12, name='epsilon')\n", "# epsilon_0 = m.Const(value=1.0006, name='epsilon_0')\n", "# mu = m.Const(value=1.2566e-6, name='mu')\n", "epsilon = 8.854e-12\n", "epsilon_0 = 1.0006\n", "mu = 1.2566e-6\n", "\n", "# Design Requirements\n", "freq_const = 13.56\n", "Do_const = 200\n", "sys_thickness = 200\n", "copper_thickness = 0.5\n", "\n", "# Define variables with bounds\n", "Freq = m.FV(value=freq_const, lb=1, ub=30,name='Freq')\n", "Do_1 = m.FV(value=Do_const, lb=100, ub=200,name='Do_1')\n", "Do_2 = m.FV(value=Do_const, lb=100, ub=200,name='Do_2')\n", "t1 = m.FV(value=copper_thickness, lb=0.5, ub=2,name='t1')\n", "t2 = m.FV(value=copper_thickness, lb=0.5, ub=2,name='t2')\n", "p1 = m.FV(value=7, lb=1, ub=100,name='p1')\n", "p2 = m.FV(value=7, lb=1, ub=100,name='p2')\n", "d = m.FV(value=3, lb=2, ub=200,name='d')\n", "N = m.FV(value=4, lb=2, ub=100, integer=True, name='N')\n", "w1 = m.FV(value=15, lb=5, ub=100,name='w1')\n", "w2 = m.FV(value=15, lb=5, ub=100,name='w2')\n", "Di_1 = m.FV(value=50, lb=5, ub=100,name='Di_1')\n", "Di_2 = m.FV(value=50, lb=5, ub=100,name='Di_2')\n", "\n", "\n", "# Variable options (0=fixed, 1=varied variable)\n", "Freq.STATUS = 0\n", "Do_1.STATUS = 0\n", "Do_2.STATUS = 0\n", "t1.STATUS = 0\n", "t2.STATUS = 0\n", "p1.STATUS = 1\n", "p2.STATUS = 1\n", "d.STATUS = 1\n", "N.STATUS = 1\n", "Di_1.STATUS = 1\n", "Di_2.STATUS = 1\n", "w1.STATUS = 1\n", "w2.STATUS = 1\n", "\n", "# Average of some variables for C_coil and L_coil equation\n", "Do = (Do_1+Do_2)/2\n", "Di = (Di_1+Di_2)/2\n", "w = (w1 + w2) / 2\n", "t = (t1 + t2) / 2\n", "\n", "# Define intermediates variables (reduce model complexity)\n", "C_coil = m.Intermediate(epsilon * epsilon_0 * np.pi * w * 1e-3 * N * (Do + Di) * 1e-3 / (2 * d * 1e-3) * (1 + d / (np.pi * w) * m.log(2 * d / (np.pi * w)) \n", " + d / (np.pi * w) * m.log(1 + 2 * d / (np.pi * w) + 2 * m.sqrt((t / d) + (t / d) ** 2))), name='C_coil')\n", "L_coil = m.Intermediate(mu * (N ** 2) * (Do + Di) * 1e-3 / 4 * m.log(2.46 * (Do + Di) / (Do - Di) + 0.2 * ((Do - Di) / (Do + Di)) ** 2), name='L_coil')\n", "resonant_freq = m.Intermediate(1/(2*np.pi*m.sqrt(L_coil*C_coil)), name='resonant_Freq')\n", "\n", "# Define constraints\n", "m.Equation(Di_1 == Di_2)\n", "# m.Equation(w1 == w2)\n", "m.Equation(Di_1 <= Do_const)\n", "m.Equation(Di_2 <= Do_const)\n", "m.Equation(t1 + t2 + d <= sys_thickness)\n", "m.Equation(p1 == (Do_1 - Di_1 - (3 * np.pi * N + 4 * np.pi - 1) / (2 * np.pi) * w1) * 2 * np.pi / (3 * np.pi * N - 1))\n", "m.Equation(p2 == (Do_2 - Di_2 - (3 * np.pi * N + 4 * np.pi - 1) / (2 * np.pi) * w2) * 2 * np.pi / (3 * np.pi * N - 1))\n", "# m.Equation(w1 >= p1)\n", "# m.Equation(w2 >= p2)\n", "m.Equation(p1 <= w1)\n", "m.Equation(p2 <= w2)\n", "m.Equation(resonant_freq <= (freq_const+0.1*freq_const)*1e6)\n", "m.Equation(resonant_freq >= (freq_const-0.1*freq_const)*1e6)\n", "\n", "# Define all variables for prediction\n", "variables = [Freq, Do_1, Do_2, d, p1, p2, N, w1, w2, Di_1, Di_2, t1, t2]\n", "\n", "# Define model from Tensorflow/Keras trained model\n", "# pred_Q_factor = Gekko_NN_TF(model, mma, m, n_output=1, activationFxn='relu').predict(variables)\n", "pred_Q_factor = m.Param(value=Gekko_NN_TF(model, mma, m, n_output=1, activationFxn='relu').predict(variables), lb=0, ub=2000, name='pred_Q_factor')\n", "\n", "# Define the GEKKO objective function\n", "m.Maximize(pred_Q_factor)\n", "\n", "m.options.SOLVER = 1 # APOPT in MINLP\n", "m.solver_options = [\n", " 'minlp_maximum_iterations 10000', \n", " 'minlp_max_iter_with_int_sol 10000', \n", " 'minlp_as_nlp 0', \n", " 'nlp_maximum_iterations 10000', \n", " 'minlp_branch_method 1', \n", " # 'minlp_integer_leaves 1',\n", " 'minlp_print_level 8', \n", " 'minlp_integer_tol 1e-4', \n", " 'minlp_gap_tol 1e-4'\n", "]\n", "m.options.IMODE = 2 # Steady-state parameter estimation mode\n", "m.solve(debug=0, disp=True)\n", "\n", "# m.open_folder() # use this to look at infeasibilities.txt\n", " \n", "# Print the results\n", "params_max_Q = []\n", "print(\"\\nOptimized design:\")\n", "for i in range(len(variables)):\n", " print(f'\\t{features[i]}: {variables[i].value[0]: .3f}')\n", " params_max_Q.append(variables[i].value)\n", "params_max_Q = np.array(params_max_Q).reshape(1, len(variables))\n", "print(f\"Maximized Q_factor: {-m.options.OBJFCNVAL}\")\n", "print(f\"Maximized Q_factor: {model.predict(params_max_Q)[0][0]: .3f}\")\n", "print('Gekko Solvetime:', m.options.SOLVETIME, 'sec')\n", "import time\n", "min_time = time.strftime(\"%H:%M:%S\", time.gmtime(m.options.SOLVETIME))\n", "print('Gekko Solvetime:', min_time)\n", "\n", "# Print L, C, Freq from the optimized design.\n", "Freq = variables[0].value[0]\n", "Do_1 = variables[1].value[0]\n", "Do_2 = variables[2].value[0]\n", "d = variables[3].value[0]\n", "p1 = variables[4].value[0]\n", "p2 = variables[5].value[0]\n", "N = variables[6].value[0]\n", "w1 = variables[7].value[0]\n", "w2 = variables[8].value[0]\n", "Di_1 = variables[9].value[0]\n", "Di_2 = variables[10].value[0]\n", "t1 = variables[11].value[0]\n", "t2 = variables[12].value[0]\n", "\n", "w = (w1+w2)/2\n", "t = (t1+t2)/2\n", "Do = (Do_1+Do_2)/2\n", "Di = (Di_1+Di_2)/2\n", "mu = 1.2566e-6\n", "epsilon = 8.854e-12\n", "epsilon_0 = 1.0006\n", "C_coil = epsilon * epsilon_0 * np.pi * w * 1e-3 * N * (Do + Di) * 1e-3 / (2 * d * 1e-3) * (1 + d / (np.pi * w) * np.log(2 * d / (np.pi * w)) + d / (np.pi * w) * np.log(1 + 2 * d / (np.pi * w) + 2 * np.sqrt((t / d) + (t / d) ** 2)))\n", "L_coil = mu * (N ** 2) * (Do + Di) * 1e-3 / 4 * np.log(2.46 * (Do + Di) / (Do - Di) + 0.2 * ((Do - Di) / (Do + Di)) ** 2)\n", "\n", "print(f'L_coil= {L_coil*10**6: .3f} uH')\n", "print(f'C_coil= {C_coil*10**12: .3f} pF')\n", "print(f'Resonant Freq = {1/(2*np.pi*np.sqrt(L_coil*C_coil))/1e6: .3f} MHz')" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.10.7" } }, "nbformat": 4, "nbformat_minor": 2 }