04 quadratic program

import%20marimo%0A%0A__generated_with%20%3D%20%220.12.0%22%0Aapp%20%3D%20marimo.App()%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20marimo%20as%20mo%0A%20%20%20%20return%20(mo%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20r%22%22%22%0A%20%20%20%20%20%20%20%20%23%20Quadratic%20program%0A%0A%20%20%20%20%20%20%20%20A%20quadratic%20program%20is%20an%20optimization%20problem%20with%20a%20quadratic%20objective%20and%0A%20%20%20%20%20%20%20%20affine%20equality%20and%20inequality%20constraints.%20A%20common%20standard%20form%20is%20the%0A%20%20%20%20%20%20%20%20following%3A%0A%0A%20%20%20%20%20%20%20%20%5C%5B%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Cbegin%7Barray%7D%7Bll%7D%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Ctext%7Bminimize%7D%20%20%20%26%20(1%2F2)x%5ETPx%20%2B%20q%5ETx%5C%5C%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Ctext%7Bsubject%20to%7D%20%26%20Gx%20%5Cleq%20h%20%5C%5C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%26%20Ax%20%3D%20b.%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Cend%7Barray%7D%0A%20%20%20%20%20%20%20%20%5C%5D%0A%0A%20%20%20%20%20%20%20%20Here%20%24P%20%5Cin%20%5Cmathcal%7BS%7D%5E%7Bn%7D_%2B%24%2C%20%24q%20%5Cin%20%5Cmathcal%7BR%7D%5En%24%2C%20%24G%20%5Cin%20%5Cmathcal%7BR%7D%5E%7Bm%20%5Ctimes%20n%7D%24%2C%20%24h%20%5Cin%20%5Cmathcal%7BR%7D%5Em%24%2C%20%24A%20%5Cin%20%5Cmathcal%7BR%7D%5E%7Bp%20%5Ctimes%20n%7D%24%2C%20and%20%24b%20%5Cin%20%5Cmathcal%7BR%7D%5Ep%24%20are%20problem%20data%20and%20%24x%20%5Cin%20%5Cmathcal%7BR%7D%5E%7Bn%7D%24%20is%20the%20optimization%20variable.%20The%20inequality%20constraint%20%24Gx%20%5Cleq%20h%24%20is%20elementwise.%0A%0A%20%20%20%20%20%20%20%20**Why%20quadratic%20programming%3F**%20Quadratic%20programs%20are%20convex%20optimization%20problems%20that%20generalize%20both%20least-squares%20and%20linear%20programming.They%20can%20be%20solved%20efficiently%20and%20reliably%2C%20even%20in%20real-time.%0A%0A%20%20%20%20%20%20%20%20**An%20example%20from%20finance.**%20A%20simple%20example%20of%20a%20quadratic%20program%20arises%20in%20finance.%20Suppose%20we%20have%20%24n%24%20different%20stocks%2C%20an%20estimate%20%24r%20%5Cin%20%5Cmathcal%7BR%7D%5En%24%20of%20the%20expected%20return%20on%20each%20stock%2C%20and%20an%20estimate%20%24%5CSigma%20%5Cin%20%5Cmathcal%7BS%7D%5E%7Bn%7D_%2B%24%20of%20the%20covariance%20of%20the%20returns.%20Then%20we%20solve%20the%20optimization%20problem%0A%0A%20%20%20%20%20%20%20%20%5C%5B%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Cbegin%7Barray%7D%7Bll%7D%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Ctext%7Bminimize%7D%20%20%20%26%20(1%2F2)x%5ET%5CSigma%20x%20-%20r%5ETx%5C%5C%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Ctext%7Bsubject%20to%7D%20%26%20x%20%5Cgeq%200%20%5C%5C%0A%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%20%26%20%5Cmathbf%7B1%7D%5ETx%20%3D%201%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20%5Cend%7Barray%7D%0A%20%20%20%20%20%20%20%20%5C%5D%0A%0A%20%20%20%20%20%20%20%20to%20find%20a%20nonnegative%20portfolio%20allocation%20%24x%20%5Cin%20%5Cmathcal%7BR%7D%5En_%2B%24%20that%20optimally%20balances%20expected%20return%20and%20variance%20of%20return.%0A%0A%20%20%20%20%20%20%20%20When%20we%20solve%20a%20quadratic%20program%2C%20in%20addition%20to%20a%20solution%20%24x%5E%5Cstar%24%2C%20we%20obtain%20a%20dual%20solution%20%24%5Clambda%5E%5Cstar%24%20corresponding%20to%20the%20inequality%20constraints.%20A%20positive%20entry%20%24%5Clambda%5E%5Cstar_i%24%20indicates%20that%20the%20constraint%20%24g_i%5ETx%20%5Cleq%20h_i%24%20holds%20with%20equality%20for%20%24x%5E%5Cstar%24%20and%20suggests%20that%20changing%20%24h_i%24%20would%20change%20the%20optimal%20value.%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20r%22%22%22%0A%20%20%20%20%20%20%20%20%23%23%20Example%0A%0A%20%20%20%20%20%20%20%20In%20this%20example%2C%20we%20use%20CVXPY%20to%20construct%20and%20solve%20a%20quadratic%20program.%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_()%3A%0A%20%20%20%20import%20cvxpy%20as%20cp%0A%20%20%20%20import%20numpy%20as%20np%0A%20%20%20%20return%20cp%2C%20np%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(%22%22%22First%20we%20generate%20synthetic%20data.%20In%20this%20problem%2C%20we%20don't%20include%20equality%20constraints%2C%20only%20inequality.%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(np)%3A%0A%20%20%20%20m%20%3D%204%0A%20%20%20%20n%20%3D%202%0A%0A%20%20%20%20np.random.seed(1)%0A%20%20%20%20q%20%3D%20np.random.randn(n)%0A%20%20%20%20G%20%3D%20np.random.randn(m%2C%20n)%0A%20%20%20%20h%20%3D%20G%20%40%20np.random.randn(n)%0A%20%20%20%20return%20G%2C%20h%2C%20m%2C%20n%2C%20q%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo%2C%20np)%3A%0A%20%20%20%20import%20wigglystuff%0A%0A%20%20%20%20P_widget%20%3D%20mo.ui.anywidget(%0A%20%20%20%20%20%20%20%20wigglystuff.Matrix(np.array(%5B%5B4.0%2C%20-1.4%5D%2C%20%5B-1.4%2C%204%5D%5D)%2C%20step%3D0.1)%0A%20%20%20%20)%0A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20f%22%22%22%0A%20%20%20%20%20%20%20%20The%20quadratic%20form%20%24P%24%20is%20equal%20to%20the%20symmetrized%20version%20of%20this%0A%20%20%20%20%20%20%20%20matrix%3A%0A%0A%20%20%20%20%20%20%20%20%7BP_widget.center()%7D%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%20P_widget%2C%20wigglystuff%0A%0A%0A%40app.cell%0Adef%20_(P_widget%2C%20np)%3A%0A%20%20%20%20P%20%3D%200.5%20*%20(np.array(P_widget.matrix)%20%2B%20np.array(P_widget.matrix).T)%0A%20%20%20%20return%20(P%2C)%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(r%22%22%22Next%2C%20we%20specify%20the%20problem.%20Notice%20that%20we%20use%20the%20%60quad_form%60%20function%20from%20CVXPY%20to%20create%20the%20quadratic%20form%20%24x%5ETPx%24.%22%22%22)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(G%2C%20P%2C%20cp%2C%20h%2C%20n%2C%20q)%3A%0A%20%20%20%20x%20%3D%20cp.Variable(n)%0A%0A%20%20%20%20problem%20%3D%20cp.Problem(%0A%20%20%20%20%20%20%20%20cp.Minimize((1%20%2F%202)%20*%20cp.quad_form(x%2C%20P)%20%2B%20q.T%20%40%20x)%2C%0A%20%20%20%20%20%20%20%20%5BG%20%40%20x%20%3C%3D%20h%5D%2C%0A%20%20%20%20)%0A%20%20%20%20_%20%3D%20problem.solve()%0A%20%20%20%20return%20problem%2C%20x%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo%2C%20problem%2C%20x)%3A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20f%22%22%22%0A%20%20%20%20%20%20%20%20The%20optimal%20value%20is%20%7Bproblem.value%3A.04f%7D.%0A%0A%20%20%20%20%20%20%20%20A%20solution%20%24x%24%20is%20%7Bmo.as_html(list(x.value))%7D%0A%20%20%20%20%20%20%20%20A%20dual%20solution%20is%20is%20%7Bmo.as_html(list(problem.constraints%5B0%5D.dual_value))%7D%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell%0Adef%20_(G%2C%20P%2C%20h%2C%20plot_contours%2C%20q%2C%20x)%3A%0A%20%20%20%20plot_contours(P%2C%20G%2C%20h%2C%20q%2C%20x.value)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(mo)%3A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20r%22%22%22%0A%20%20%20%20%20%20%20%20In%20this%20plot%2C%20the%20gray%20shaded%20region%20is%20the%20feasible%20region%20(points%20satisfying%20the%20inequality)%2C%20and%20the%20ellipses%20are%20level%20curves%20of%20the%20quadratic%20form.%0A%0A%20%20%20%20%20%20%20%20**%F0%9F%8C%8A%20Try%20it!**%20Try%20changing%20the%20entries%20of%20%24P%24%20above%20with%20your%20mouse.%20How%20do%20the%0A%20%20%20%20%20%20%20%20level%20curves%20and%20the%20optimal%20value%20of%20%24x%24%20change%3F%20Can%20you%20explain%20what%20you%20see%3F%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(P%2C%20mo)%3A%0A%20%20%20%20mo.md(%0A%20%20%20%20%20%20%20%20rf%22%22%22%0A%20%20%20%20%20%20%20%20The%20above%20contour%20lines%20were%20generated%20with%0A%0A%20%20%20%20%20%20%20%20%5C%5B%0A%20%20%20%20%20%20%20%20P%3D%20%5Cbegin%7B%7Bbmatrix%7D%7D%0A%20%20%20%20%20%20%20%20%7BP%5B0%2C%200%5D%3A.01f%7D%20%26%20%7BP%5B0%2C%201%5D%3A.01f%7D%20%5C%5C%0A%20%20%20%20%20%20%20%20%7BP%5B1%2C%200%5D%3A.01f%7D%20%26%20%7BP%5B1%2C%201%5D%3A.01f%7D%20%5C%5C%0A%20%20%20%20%20%20%20%20%5Cend%7B%7Bbmatrix%7D%7D%0A%20%20%20%20%20%20%20%20%5C%5D%0A%20%20%20%20%20%20%20%20%22%22%22%0A%20%20%20%20)%0A%20%20%20%20return%0A%0A%0A%40app.cell(hide_code%3DTrue)%0Adef%20_(np)%3A%0A%20%20%20%20def%20plot_contours(P%2C%20G%2C%20h%2C%20q%2C%20x_star)%3A%0A%20%20%20%20%20%20%20%20import%20matplotlib.pyplot%20as%20plt%0A%0A%20%20%20%20%20%20%20%20%23%20Create%20a%20grid%20of%20x%20and%20y%20values.%0A%20%20%20%20%20%20%20%20x%20%3D%20np.linspace(-5%2C%205%2C%20400)%0A%20%20%20%20%20%20%20%20y%20%3D%20np.linspace(-5%2C%205%2C%20400)%0A%20%20%20%20%20%20%20%20X%2C%20Y%20%3D%20np.meshgrid(x%2C%20y)%0A%0A%20%20%20%20%20%20%20%20%23%20Compute%20the%20quadratic%20form%20Q(x%2C%20y)%20%3D%20a*x%5E2%20%2B%202*b*x*y%20%2B%20c*y%5E2.%0A%20%20%20%20%20%20%20%20%23%20Here%2C%20a%20%3D%20P%5B0%2C0%5D%2C%20b%20%3D%20P%5B0%2C1%5D%20(and%20P%5B1%2C0%5D)%2C%20c%20%3D%20P%5B1%2C1%5D%0A%20%20%20%20%20%20%20%20Z%20%3D%20(%0A%20%20%20%20%20%20%20%20%20%20%20%200.5%20*%20(P%5B0%2C%200%5D%20*%20X**2%20%2B%202%20*%20P%5B0%2C%201%5D%20*%20X%20*%20Y%20%2B%20P%5B1%2C%201%5D%20*%20Y**2)%0A%20%20%20%20%20%20%20%20%20%20%20%20%2B%20q%5B0%5D%20*%20X%0A%20%20%20%20%20%20%20%20%20%20%20%20%2B%20q%5B1%5D%20*%20Y%0A%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20%23%20---%20Evaluate%20the%20constraints%20on%20the%20grid%20---%0A%20%20%20%20%20%20%20%20%23%20We%20stack%20X%20and%20Y%20to%20get%20a%20list%20of%20(x%2Cy)%20points.%0A%20%20%20%20%20%20%20%20points%20%3D%20np.vstack(%5BX.ravel()%2C%20Y.ravel()%5D).T%0A%0A%20%20%20%20%20%20%20%20%23%20Start%20with%20all%20points%20feasible%0A%20%20%20%20%20%20%20%20feasible%20%3D%20np.ones(points.shape%5B0%5D%2C%20dtype%3Dbool)%0A%0A%20%20%20%20%20%20%20%20%23%20Apply%20the%20inequality%20constraints%20Gx%20%3C%3D%20h.%0A%20%20%20%20%20%20%20%20%23%20Each%20row%20of%20G%20and%20corresponding%20h%20defines%20a%20condition.%0A%20%20%20%20%20%20%20%20for%20i%20in%20range(G.shape%5B0%5D)%3A%0A%20%20%20%20%20%20%20%20%20%20%20%20%23%20For%20a%20given%20point%20x%2C%20the%20condition%20is%3A%20G%5Bi%2C0%5D*x%20%2B%20G%5Bi%2C1%5D*y%20%3C%3D%20h%5Bi%5D%0A%20%20%20%20%20%20%20%20%20%20%20%20feasible%20%26%3D%20points.dot(G%5Bi%5D)%20%3C%3D%20h%5Bi%5D%20%2B%201e-8%20%20%23%20small%20tolerance%0A%20%20%20%20%20%20%20%20%23%20Reshape%20the%20boolean%20mask%20back%20to%20grid%20shape.%0A%20%20%20%20%20%20%20%20feasible_grid%20%3D%20feasible.reshape(X.shape)%0A%0A%20%20%20%20%20%20%20%20%23%20---%20Plot%20the%20feasible%20region%20and%20contour%20lines---%0A%20%20%20%20%20%20%20%20plt.figure(figsize%3D(8%2C%206))%0A%0A%20%20%20%20%20%20%20%20%23%20Use%20contourf%20to%20fill%20the%20region%20where%20feasible_grid%20is%20True.%0A%20%20%20%20%20%20%20%20%23%20We%20define%20two%20levels%2C%20so%20that%20points%20that%20are%20True%20(feasible)%20get%20one%0A%20%20%20%20%20%20%20%20%23%20color.%0A%20%20%20%20%20%20%20%20plt.contourf(%0A%20%20%20%20%20%20%20%20%20%20%20%20X%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20Y%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20feasible_grid%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20levels%3D%5B-0.5%2C%200.5%2C%201.5%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20colors%3D%5B%22white%22%2C%20%22gray%22%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20alpha%3D0.5%2C%0A%20%20%20%20%20%20%20%20)%0A%0A%20%20%20%20%20%20%20%20contours%20%3D%20plt.contour(X%2C%20Y%2C%20Z%2C%20levels%3D10%2C%20cmap%3D%22viridis%22)%0A%20%20%20%20%20%20%20%20plt.clabel(contours%2C%20inline%3DTrue%2C%20fontsize%3D8)%0A%20%20%20%20%20%20%20%20plt.title(%22Feasible%20region%20and%20level%20curves%22)%0A%20%20%20%20%20%20%20%20plt.xlabel(%22%24x_1%24%22)%0A%20%20%20%20%20%20%20%20plt.ylabel(%22%24y_2%24%22)%0A%20%20%20%20%20%20%20%20%23%20plt.colorbar(contours%2C%20label%3D'Q(x%2C%20y)')%0A%0A%20%20%20%20%20%20%20%20ax%20%3D%20plt.gca()%0A%20%20%20%20%20%20%20%20%23%20Optionally%2C%20mark%20and%20label%20the%20point%20x_star.%0A%20%20%20%20%20%20%20%20ax.plot(x_star%5B0%5D%2C%20x_star%5B1%5D%2C%20%22ko%22%2C%20markersize%3D5)%0A%20%20%20%20%20%20%20%20ax.text(%0A%20%20%20%20%20%20%20%20%20%20%20%20x_star%5B0%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20x_star%5B1%5D%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20r%22%24%5Cmathbf%7Bx%7D%5E%5Cstar%24%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20color%3D%22black%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20fontsize%3D12%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20verticalalignment%3D%22bottom%22%2C%0A%20%20%20%20%20%20%20%20%20%20%20%20horizontalalignment%3D%22right%22%2C%0A%20%20%20%20%20%20%20%20)%0A%20%20%20%20%20%20%20%20return%20plt.gca()%0A%20%20%20%20return%20(plot_contours%2C)%0A%0A%0Aif%20__name__%20%3D%3D%20%22__main__%22%3A%0A%20%20%20%20app.run()%0A