Ask Your Question
1

Creating a function in Sage

asked 2020-10-14 17:29:52 +0100

florin gravatar image

Hi I have a sequence of commands which does Lagrangian optimization

var('a, b, x, y, lam')
assume(a>=0,b>=0,x>= 0)
f = x*y #minimiser et maximiser sous contrainte d'egalite h=0
h = x^2/a^2 + y^2/b^2-1
L = f + h * lam
DL = L.gradient([x, y])
sol = solve([DL[k] == 0 for k in range(len(DL))]+[ h == 0],x, y, lam,
solution_dict=True)
ti=[["$\\lambda$","$x$", "$y$", "$f(x, y)$"]]
re=[(sol[j]
[lam],sol[j][x], sol[j][y],f(x=sol[j][x], y=sol[j][y])) for j in
range(4)]
table(ti + re, header_row = True)

and I would like to turn it into a function. I add the declaration, but I get no output. I know neither why, nor how to debug this. Thanks for any help :)

def optL(f,h):
       var('a, b, x, y, lam')
       assume(a>=0,b>=0,x>= 0)
       L = f + h * lam
       DL = L.gradient([x, y])
       sol = solve([DL[k] == 0 for k in range(len(DL))]+[ h == 0],x, y, lam,
       solution_dict=True)
       ti=[["$\\lambda$","$x$", "$y$", "$f(x, y)$"]]
       re=[(sol[j][lam],sol[j][x], sol[j][y],f(x=sol[j][x], y=sol[j][y])) for j in range(4)]
       table(ti + re, header_row = True)
optL(x*y,x^2/a^2 + y^2/b^2-1)
edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
0

answered 2020-10-14 17:57:39 +0100

Emmanuel Charpentier gravatar image

updated 2020-10-14 17:58:59 +0100

Without looking at what your code does, a couple remarks :

  1. You need to return the value(s) of the function

  2. Returning atable isn't probably what you want ; what about a dictionary ?

edit flag offensive delete link more

Comments

Merci, Emmanuel, now it works :)

def optL(f,h):
   var('a, b, x, y, lam')
   assume(a>=0,b>=0,x>= 0)
   L = f + h * lam
   DL = L.gradient([x, y])
   sol = solve([DL[k] == 0 for k in range(len(DL))]+[ h == 0],x, y, lam,
   solution_dict=True)
   ti=[["$\\lambda$","$x$", "$y$", "$f(x, y)$"]]
   re=[(sol[j][lam],sol[j][x], sol[j][y],f(x=sol[j][x], y=sol[j][y])) for j in range(len(sol))]
   ans=table(ti + re, header_row = True)
   show(ans)
   return re
ex5=optL(x*y,x^2/a^2+ y^2/b^2-1)
ex5
florin gravatar imageflorin ( 2020-10-14 20:07:50 +0100 )edit

I turned @Emmanuel Charpentier's comment into an answer so you can accept it and thus mark the question as solved.

slelievre gravatar imageslelievre ( 2021-01-01 21:49:10 +0100 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2020-10-14 17:29:52 +0100

Seen: 620 times

Last updated: Oct 14 '20