From d1840a9ed7abc7f96570092158c14829daf44792 Mon Sep 17 00:00:00 2001 From: youainti Date: Thu, 23 Sep 2021 11:13:01 -0700 Subject: [PATCH] Got a basic optimality condition function working --- Code/README.md | 6 +- Code/connect_transition_to_optimality.ipynb | 179 ++++++++++++++------ 2 files changed, 133 insertions(+), 52 deletions(-) diff --git a/Code/README.md b/Code/README.md index 2eff7f9..7baeb4a 100644 --- a/Code/README.md +++ b/Code/README.md @@ -1,9 +1,11 @@ # COMPUTATIONAL TODO -## Next steps +## Completed steps - implement 'launch function as a function' portion - - This will probably use a neural network as a function type approach. - substitute the transition functions into the optimality conditions. + +## Next steps +- create the iterated optimality conditions - use these optimality conditions to create a loss function - add boundary conditions to loss function - get a basic gradient descent/optimization of launch function working. diff --git a/Code/connect_transition_to_optimality.ipynb b/Code/connect_transition_to_optimality.ipynb index 015dc22..cb70164 100644 --- a/Code/connect_transition_to_optimality.ipynb +++ b/Code/connect_transition_to_optimality.ipynb @@ -3,7 +3,7 @@ { "cell_type": "code", "execution_count": 1, - "id": "behavioral-session", + "id": "expired-austria", "metadata": { "tags": [] }, @@ -16,7 +16,7 @@ }, { "cell_type": "markdown", - "id": "vocational-rover", + "id": "damaged-accountability", "metadata": {}, "source": [ "# Setup Functions\n", @@ -26,7 +26,7 @@ { "cell_type": "code", "execution_count": 2, - "id": "physical-guidance", + "id": "modular-memorabilia", "metadata": {}, "outputs": [], "source": [ @@ -73,7 +73,7 @@ }, { "cell_type": "markdown", - "id": "difficult-drinking", + "id": "direct-picture", "metadata": {}, "source": [ "# functions related to transitions" @@ -81,12 +81,12 @@ }, { "cell_type": "code", - "execution_count": 3, - "id": "beautiful-northwest", + "execution_count": 37, + "id": "bridal-ordinary", "metadata": {}, "outputs": [], "source": [ - "def single_transition(item_to_iterate, laws_motion, profit, stocks, debris ):\n", + "def single_transition(item_to_iterate, laws_motion_fn, profit_fn, stocks, debris, launch_fn):\n", " \"\"\"\n", " This function represents the inverted envelope conditions.\n", " It allows us to describe the derivatives of the value function evaluated at time $t+1$ in terms based in time period $t$.\n", @@ -101,20 +101,27 @@ " #it consists of the derivative of the laws of motion with respect to stocks and debris\n", " \n", " #Get the jacobian\n", - " a = jacobian(laws_motion, (stocks,debris))\n", + " a = jacobian(laws_motion_fn, (stocks,debris, launch_fn(stocks,debris)))\n", " \n", " #Reassemble the Jacobian nested tuples into the appropriate tensor\n", " A = BETA * torch.cat((torch.cat((a[0][0],a[0][1]),dim=1),torch.cat((a[1][0],a[1][1]),dim=1)), dim=0)\n", - " \n", + "\n", " #TODO: figure out some diagnostics for this section\n", " #Possibilities include:\n", " # - Det(A) ~= 0\n", " # - EigVal(A) ~= 0\n", " # - A.inverse() with a try catch system to record types of returns\n", " #Alternatively, \n", + " #if abs(a.det())\n", " \n", " #Calculate the item to transition\n", - " T = item_to_iterate - torch.cat(jacobian(profit,(stocks, debris))) \n", + " f_jacobians = jacobian(profit_fn,(stocks, debris, launch_fn(stocks,debris)))\n", + "\n", + " #issue with shape here: my launch function is for all launches, not just a single launch.\n", + " f_theta = torch.cat([f_jacobians[0][0], f_jacobians[1][0]],axis=0) \n", + "\n", + " T = item_to_iterate - f_theta\n", + "\n", " #Includes rearranging the jacobian of profit.\n", "\n", " #Return the transitioned values\n", @@ -125,30 +132,31 @@ " \"\"\"\n", " \"\"\"\n", " #unpack states and functions\n", - " stocks, debris,profit, laws_motion, item_to_transition = data_in\n", + " stocks, debris,profit_fn, laws_motion_fn, item_to_transition,launch_fn = data_in\n", " \n", " #Calculate new states\n", - " new_stocks, new_debris = laws_motion(stocks,debris)\n", + " new_stocks, new_debris = laws_motion_fn(stocks,debris, launch_fn(stocks,debris))\n", " \n", " #WARNING: RECURSION: You may break your head...\n", " #This gets the transition of the value function derivatives over time.\n", " transitioned = single_transition(\n", " item_to_transition, #item to iterate, i.e. the derivatives of the value function\n", " #functions\n", - " laws_motion, \n", - " profit, \n", - " stocks, debris #states\n", + " laws_motion_fn, \n", + " profit_fn, \n", + " stocks, debris, #states\n", + " launch_fn #launch function\n", " )\n", " \n", " #collects the data back together for return, including the updated state variables\n", - " data_out = new_stocks, new_debris, profit, laws_motion, transitioned\n", + " data_out = new_stocks, new_debris, profit_fn, laws_motion_fn, transitioned, launch_fn\n", " \n", " return data_out" ] }, { "cell_type": "markdown", - "id": "entertaining-theorem", + "id": "miniature-karaoke", "metadata": {}, "source": [ "## Setup functions related to the problem" @@ -156,24 +164,13 @@ }, { "cell_type": "code", - "execution_count": 4, - "id": "modern-kentucky", - "metadata": {}, + "execution_count": 38, + "id": "bright-minimum", + "metadata": { + "tags": [] + }, "outputs": [], "source": [ - "### Classes\n", - "class States():\n", - " \"\"\"\n", - " This class represents the state variables associated with the problems.\n", - " \n", - " In this problem, the two types of states are constellation stocks and debris.\n", - " \n", - " I'm not sure how useful it will be. We'll see. It is missing a lot\n", - " \"\"\"\n", - " def __init__(self, satellite_stock, debris):\n", - " self.stock = satellite_stock\n", - " self.debris = debris\n", - " \n", " \n", "\n", "### functions\n", @@ -192,34 +189,33 @@ " \"\"\"\n", " return torch.ones(5, requires_grad=True)\n", "\n", - "def laws_of_motion(stock, debris):\n", + "def laws_of_motion(stock, debris, launch):\n", " \"\"\"\n", " This function updates state variables (stock and debris), according \n", " to the laws of motion.\n", " \n", " It returns the state variables as \n", " \"\"\"\n", - " l = launches(stock,debris)\n", - " #Notes: Launches is a global function.\n", + "\n", " s = survival(stock,debris)\n", " #Notes: Survival is a global function.\n", " \n", - " new_stock = stock*s + l\n", + " new_stock = stock*s + launch\n", " \n", " \n", " #TODO: Currently Ignoring autocatalysis\n", - " new_debris = (1-DELTA)*debris + LAUNCH_DEBRIS_RATE * l.sum() + COLLISION_DEBRIS_RATE*(1-s) @ stock\n", + " new_debris = (1-DELTA)*debris + LAUNCH_DEBRIS_RATE * launch.sum() + COLLISION_DEBRIS_RATE*(1-s) @ stock\n", " \n", " return (new_stock, new_debris)\n", "\n", "#This is not a good specification of the profit function, but it will work for now.\n", - "def profit(stock, debris):\n", - " return UTIL_WEIGHTS @ stock" + "def profit(stock, debris, launches):\n", + " return UTIL_WEIGHTS @ stock - LAUNCH_COST*launches" ] }, { "cell_type": "markdown", - "id": "broadband-technique", + "id": "conservative-ukraine", "metadata": {}, "source": [ "# Actual calculations" @@ -227,8 +223,8 @@ }, { "cell_type": "code", - "execution_count": 5, - "id": "adjustable-harvey", + "execution_count": 39, + "id": "initial-mathematics", "metadata": {}, "outputs": [], "source": [ @@ -245,13 +241,15 @@ "#CHANGE LATER: Launch is currently a value, should be a function (i.e. neural network)\n", "launches = test_launch\n", "\n", - "#compose the functions together.\n", - "base_data = (stocks,debris, profit, laws_of_motion, torch.ones(6, requires_grad=True))\n", + "#Starting point\n", + "# Stocks, debris, profit fn, laws of motion, \n", + "base_data = (stocks,debris, profit, laws_of_motion, torch.ones(6, requires_grad=True),launches)\n", "\n", "#Parameters\n", "SCALING = torch.ones(5)\n", - "DELTA = 0.9 \n", + "DELTA = 0.9\n", "LAUNCH_DEBRIS_RATE = 0.005\n", + "LAUNCH_COST = 1.0\n", "COLLISION_DEBRIS_RATE = 0.0007\n", "UTIL_WEIGHTS = torch.tensor([1,-0.2,0,0,0])\n", "BETA = 0.95" @@ -259,8 +257,8 @@ }, { "cell_type": "code", - "execution_count": 6, - "id": "cordless-wages", + "execution_count": 75, + "id": "nuclear-definition", "metadata": {}, "outputs": [ { @@ -304,16 +302,97 @@ }, { "cell_type": "markdown", - "id": "shaped-zambia", + "id": "casual-annex", "metadata": {}, "source": [ "Also, maybe I can create a `Model` class that upon construction will capture the necesary constants, functions, etc.\n" ] }, + { + "cell_type": "markdown", + "id": "stopped-socket", + "metadata": {}, + "source": [ + "# Optimatility conditions" + ] + }, + { + "cell_type": "code", + "execution_count": 193, + "id": "excessive-script", + "metadata": {}, + "outputs": [], + "source": [ + "#Optimality condition\n", + "def optimality(stocks,debris,profit_fn,laws_motion_fn,launch_fn, iterated_item):\n", + " #Derivative of the value function with respect to choice functions\n", + " #this returns derivatives with respect to every launch, so I've removed that\n", + " fx = jacobian(profit_fn, (stocks,debris,launch_fn(stocks,debris)))[-1][:,0]\n", + " \n", + " \n", + " #The following returns a tuple of tuples of tensors.\n", + " #the first tuple contains jacobians related to laws of motion for stocks\n", + " #the second tuple contains jacobians related to laws of motion for debris.\n", + " #we need the derivatives related to both\n", + " b = jacobian(laws_of_motion,(stocks,debris,launches(stocks,debris)))\n", + " B = torch.cat((b[0][2],b[1][2].T),axis=1)\n", + "\n", + "\n", + " return fx + BETA * B @ iterated_item" + ] + }, + { + "cell_type": "code", + "execution_count": 195, + "id": "unlikely-coverage", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "tensor(49.4968)" + ] + }, + "execution_count": 195, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "sum(optimality(stocks,debris,profit,laws_of_motion,launches,tmp_result)**2)" + ] + }, + { + "cell_type": "markdown", + "id": "endless-occupation", + "metadata": {}, + "source": [ + "## Now to set up the recursive set of optimatliy conditions" + ] + }, + { + "cell_type": "code", + "execution_count": 179, + "id": "valuable-bleeding", + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "tensor([0.0300, 2.0300, 3.0300, 4.0300, 5.0300])" + ] + }, + "execution_count": 179, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [] + }, { "cell_type": "code", "execution_count": null, - "id": "canadian-excitement", + "id": "subjective-chassis", "metadata": {}, "outputs": [], "source": []