{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "## CS 132 Homework 5 A \n", "\n", "### Due Thursday August 12th at Midnight (1 minute after 11:59pm) in Gradescope (with grace period of 6 hours)\n", "### Homeworks may be submitted up to 24 hours late with a 10% penalty (same grace period)\n", "\n", "\n", "Please read through the entire notebook, reading the expository material and then do the problems, following the instuctions to try various features; among the exposition of various features of Python there\n", "are three problems which will be graded. All problems are worth 10 points. \n", "\n", "You will need to complete this notebook and then convert to a PDF file in order to submit to Gradescope. Instructions are given here:\n", "\n", "https://www.cs.bu.edu/fac/snyder/cs132/HWSubmissionInstructions.html\n", "\n" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# Here are some imports which will be used in the code in the rest of the lab \n", "\n", "# Imports used for the code in CS 237\n", "\n", "import numpy as np # arrays and functions which operate on arrays\n", "\n", "import matplotlib.pyplot as plt # normal plotting\n", "\n", "\n", "# NOTE: You may not use any other libraries than those listed here without explicit permission." ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "# Gaussian Elimination\n", "\n", "# number of digits of precision to print out\n", "\n", "prec = 4\n", "\n", "########################################################################################################\n", "\n", "# Returns the Row Echelon (upper triangular) form\n", "\n", "def forwardElimination(A,traceLevel=0):\n", " \n", " A = (np.copy(A)).astype(float)\n", " \n", " if (traceLevel > 0):\n", " print(\"Running Forward Elimination on:\\n\")\n", " print(np.round(A, decimals=4),\"\\n\")\n", " print()\n", " \n", " (numRows,numCols) = A.shape\n", " \n", " # r = row we are currently working on pivot value is A[r][c]\n", " r = 0 \n", " for c in range(numCols): # solve for variable in column c \n", " # find row in column c with first non-zero element, and exchange with row r \n", " r1 = r\n", " while(r1 < numRows):\n", " if (not np.isclose(A[r1][c],0.0)): # A[r1][c] is first non-zero element at or below r in column c\n", " break\n", " r1 += 1\n", " \n", " if(r1 == numRows): # all zeros below r in this column\n", " continue # go on to the next column, but still working on same row \n", " \n", " if(r1 != r): \n", " # exchange rows r1 and r\n", " A[[r1,r],:] = A[[r,r1],:] \n", " if (traceLevel == 2): \n", " print(\"Exchange R\" + str(r1+1) + \" and R\" + str(r+1) + \"\\n\") \n", " print(np.round(A, decimals=4)) \n", " print() \n", "\n", " # now use pivot A[r][c] to eliminate all vars in this column below row r\n", " for r2 in range(r+1,numRows):\n", " rowReduce(A,r,c,r2,traceLevel)\n", " \n", " r += 1 \n", " if (r >= numRows):\n", " break\n", " \n", " return A\n", "\n", "# for pivot A[r][c], eliminate variable in location A[r2][c] in row r2 using row operations\n", "\n", "def rowReduce(A,r,c,r2,traceLevel=0):\n", "\n", " if(not np.isclose(A[r2][c],0.0)):\n", "\n", " factor = -A[r2][c]/A[r][c] \n", " A[r2] += factor * A[r]\n", " \n", " if(traceLevel == 2):\n", " print(\"R\" + str(r2+1) + \" += \" + str(np.around(factor,prec)) + \"*R\" + str(r+1) + \"\\n\") \n", " print(np.round(A, decimals=4))\n", " print()\n", "\n", "# Take a Row Echelon Form and return a Reduced Row Echelon Form\n", "\n", "def backwardSubstitute(A,augmented=True,traceLevel=0): \n", " \n", " numRows,numCols = A.shape\n", " \n", " if (A.dtype != 'float64'):\n", " A = A.astype(float)\n", "\n", " # now back-substitute the variables from bottom row to top\n", " if (traceLevel > 0):\n", " print(\"Creating Reduced Row Echelon Form...\\n\") \n", "\n", " for r in range(numRows):\n", "\n", " # find variable in this row\n", " for c in range(numCols):\n", " if(not np.isclose(A[r][c],0.0)):\n", " break \n", " \n", " if ((augmented and c >= numCols-1) or (c >= numCols)): # inconsistent or redundant row\n", " continue \n", " \n", " # A[r][c] is variable to eliminate\n", " \n", " factor = A[r][c]\n", " \n", " if (np.isclose(factor,0.0)): # already eliminated\n", " continue\n", " \n", " if(not np.isclose(factor,1.0)): \n", " A[r] *= 1/factor\n", " if (traceLevel == 2):\n", " print(\"R\" + str(r+1) + \" = R\" + str(r+1) + \"/\" + str(np.around(factor,prec)) + \"\\n\") \n", " print(np.round(A, decimals=4))\n", " print()\n", "\n", " for r2 in range(r): \n", " rowReduce(A,r,c,r2,traceLevel)\n", " \n", " return A \n", "\n", " \n", "# try to find row of all zeros except for last column, in augmented matrix. \n", "\n", "def noSolution(A):\n", " numRows,numCols = A.shape\n", " for r in range(numRows-1,-1,-1): # start from bottom, since inconsistent rows end up there\n", " for c in range(numCols):\n", " if(not np.isclose(A[r][c],0.0)): # found first non-0 in this row\n", " if(c == numCols-1):\n", " return True\n", " else:\n", " break\n", " return False\n", "\n", "########################################################################################################\n", "\n", "# Runs GE and returns a reduced row echelon form\n", "\n", "# If b == None assumes that A is NOT an augmented matrix, and runs GE and returns Reduced Row Echelon Form\n", "\n", "# If b is a column matrix then adjoins it to A and runs GE;\n", "# Always returns the Reduced Row Echelon Form\n", "# If inconsistent then also prints out \"Inconsistent!\"\n", "\n", "# If b is a length n array instead of a 1 x n array (column vector)\n", "# b will be converted to a column vector, for convenience. \n", "\n", "# traceLevel 0 will not print anything out during the run\n", "# traceLevel 1 will print out various stages of the process and the intermediate matrices\n", "# traceLevel 2 will also print out the row operations used at each step. \n", "\n", "# If you want to produce an Echelon Form (NOT reduced), then use forwardElimination instead. \n", "\n", "# See examples for more explanation of how to use this\n", "\n", "def GaussianElimination(A,b=None, traceLevel = 0):\n", " if( type(b) != type(None)):\n", " if( (A.shape)[0] == 1 ):\n", " b = np.array([b]).T\n", " Ab = (np.hstack((A.copy(),b))).astype(float)\n", " else:\n", " Ab = A.copy().astype(float)\n", " \n", " if( traceLevel > 0 and type(b) == type(None)):\n", " print(\"Creating Reduced Row Echelon Form:\\n\")\n", " print(np.round(Ab, decimals=4))\n", " print()\n", " elif( traceLevel > 0 ):\n", " print(\"Running Gaussian Elimination on augmented matrix:\\n\")\n", " print(np.round(Ab, decimals=4))\n", " print()\n", " \n", " B = forwardElimination(Ab,traceLevel)\n", " \n", " if( traceLevel > 0 ):\n", " print(\"Echelon Form:\\n\")\n", " print( np.round(B, decimals=4) + np.zeros(B.shape),\"\\n\") # adding -0.0 + 0 gives 0.0\n", " print()\n", " \n", " if ( type(b) != type(None) and noSolution(B) ):\n", " print(\"No solution!\")\n", "\n", " C = backwardSubstitute(B,type(b)!=type(None),traceLevel)\n", " \n", " if( traceLevel > 0 ):\n", " print(\"Reduced Row Echelon Form:\\n\")\n", " print( np.round(C, decimals=4) + np.zeros(C.shape),\"\\n\") # adding -0.0 + 0 gives 0.0\n", " print()\n", " \n", " return C\n", "\n", "########################################################################################################\n", " \n", "def GE(A,b=None, traceLevel = 0):\n", " return GaussianElimination(A,b,traceLevel)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Problem One\n", "\n", "\n", "\n", "A kindergartner is going through that phase when she only eats certain\n", "foods. For breakfast she only eats Lucky Charms, Cap'n Crunch, or\n", "Froot Loops. Her parent notices that if she eats one cereal on a given\n", "day, on the next day she will eat the same cereal with probability of\n", "50\\%, and she will choose one of the other cereals with probability of\n", "25\\%.\n", "\n", "(A) What is the stochastic matrix for this situation? \n", "\n", "(B) If she chooses Cap'n Crunch on Monday, what is the\n", " probability she will choose Froot Loops on Wednesday?\n", " \n", "(C) Consider this question: \"If she chooses Lucky Charms on Sunday and Cap'n Crunch on Monday, what is the\n", " probability she will choose Froot Loops on Wednesday?\" The answer to this question is the same as the answer to (B). Explain carefully why this is the case. \n", "\n", "(D) If the parent would like to predict what cereal the child will eat on some day in the\n", "very distant future, what are the probabilities for each cereal?\n", "\n", "For (D), you must run the experiment until it reaches the steady-state vector. Show a trace of the experiment (do not run the experiment much past the steady-state solution, to save space. \n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution:**\n", " " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Problem Two\n", "\n", "The weather in Boston is either good, indifferent, or\n", " bad on any given day. If the weather is good today, there is a 60\\%\n", " chance the weather will be good tomorrow, a 30\\% chance the weather\n", " will be indifferent, and a 10\\% chance the weather will be bad. If\n", " the weather is indifferent today, it will be good tomorrow with\n", " probability 0.40 and indifferent with probability 0.30. Finally, if\n", " the weather is bad today it will be good tomorrow with probability\n", " 0.40 and indifferent with probability 0.50. \n", " \n", "(A) Give the stochastic matrix for this problem. Order the\n", " weather types as: (good, indifferent, bad).\n", " \n", "(B) Suppose there is a 50\\% chance of good weather today and a 50\\%\n", " chance of indifferent weather. What are the chances of bad weather\n", " tomorrow?\n", " \n", "(C) What is the long-term forecast for weather in Boston? Solve this by using\n", "the algebraic technique shown in lecture for calculating the steady-state\n", "vector directly from the stochastic matrix. \n", " \n", "(D) Run an experiment in Python to determined the steady-state vector; it should, of course, be very\n", "close to the result found in (C). Show the trace of the experiment. " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "**Solution:**\n", "\n" ] } ], "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.7.6" } }, "nbformat": 4, "nbformat_minor": 2 }