Project 1.5a_optimization
system:sage


<h1>Part I: An interactive program for Unconstrained Optimization</h1>
<p><br /> The second derivative test generalizes naturally to give a test for multivariable optimization. If $Df$ represents the matrix derivative of a function $f(x,y)$, then critical points are where $Df(x,y)=0$ . The second derivative, $D2f(x,y)$, called the Hessian matrix, can be used to test if a critical point is a maximum (both eigenvalues negative), minimum (both eigenvalues positive), or a saddle point (eigenvalues have different signs). The eigenvectors tell the directions in which the second derivative is a maximum or minimum. <br /><br /> In this project you are supposed to <br /> a) Extend this sage notebook so the user can specify any two dimensional objective function <br /> b) Develop the same application in MatLab using a GUI interface <br /> c) Develop the same application in python using a GUI interface</p>

{{{id=3|
def illustrate_test(f, fvals=(-5,5), tvals=(-3,3)):
    """
    Return a table of critical points, eigenvalues, and eigenvectors, plus a 3D graph illustrating the curvature at each critical point.

    INPUTS:
      f -- a function.  Usually this should be a callable function, defined like f(x,y)=x^2+y^2, for example.
      fvals (default (-5,5)) -- minimum and maximum values for the function plot
      tvals (default (-3,3)) -- minimum and maximum values for the parameter in the curves at each critical point
    """
    line_thickness=3
    surface_color='blue'
    line_color = ['red','blue']
    x,y,t=var('x,y,t')
    tmin, tmax = tvals
    fmin, fmax = fvals
    f_parametric = vector([x,y,f(x,y)])
    picture = plot3d(f,(x,fmin,fmax),(y,fmin,fmax), opacity=0.3, surface_color=surface_color)
    hessian=f(x,y).hessian()
    solns=solve([eqn(x,y)==0 for eqn in f.gradient()],[x,y], solution_dict=True)

    table=[['Color','Critical Point', 'Eigenvalue', 'Eigenvector']]
    pictures=[]
    i=0
    for solution in [s for s in solns if s[x] in RR and s[y] in RR]:
        soln = dict([(key,RDF(val)) for key,val in solution.items()])
        critical_point = f_parametric(soln)
        hessian_at_point = hessian(soln).change_ring(RDF)
        evals,evecs=hessian_at_point.eigenmatrix_right()
        evals = evals.diagonal()
        picture += point3d(critical_point.list(), color='green', size=10)
        for eval,evec in zip(evals,evecs.columns()):
            direction_plot = evec.plot(width=line_thickness,color=line_color[i%2]).plot3d().translate(critical_point)    
            curve_formula = f_parametric(x=soln[x]+t*evec[0],y=soln[y]+t*evec[1])
            curve = parametric_plot(curve_formula, (t,tmin,tmax),color=line_color[i%2],thickness=line_thickness)
            picture+=curve+direction_plot
            pictures += [curve]
            table.append(["<span style='color: %(c)s'>%(c)s</span>"%{'c':line_color[i%2]},critical_point, eval, (evec)])
            i+=1

    return table, picture

var('x,y')
@interact
def _(f=selector([-x^2+y^2, x*y*exp((-x^2-y^2)/3), x^3-3*x+y^2-4*y], label="$f(x,y)=$")):
    f(x,y)=f
    table, graph=illustrate_test(f)
    html.table(table,header=True)
    show(graph)
///
}}}

{{{id=4|

///
}}}