(%i11)
Gauss(quad,variable):=block([q,var,liste,i,a,l], q:quad,var:variable,liste:[],display(q,var,liste), while var#[] do ( a:diff(q,var[1],2)/2,l:subst(0,var[1],diff(q,var[1])),display(a,l), if a=0 and l=0 then var:rest(var) elseif a#0 then ( liste:append(liste,[[a,var[1]+l/(2*a)]]),q:subst(0,var[1],q)-l^2/(4*a),var:rest(var)) else ( i:2,while subst(0,var[i],diff(l,var[i]))=0 do i:i+1,display(var[i]), a:diff(q,var[i],2)/2,display(a), if a#0 then ( liste:append(liste,[[a,var[i]+l/(2*a)]]),q:subst(0,var[i],q)-l^2/(4*a),var:append(rest(var,i-1-length(var)),rest(var,i))) else ( a:diff(q,var[1],1,var[i],1),display(a), l1:subst([var[1]=0,var[i]=0],diff(q,var[1])),li:subst([var[1]=0,var[i]=0],diff(q,var[i])),display(l1,li), liste:append(liste,[[a/4,var[1]+var[i]+(l1+li)/a],[-a/4,var[1]-var[i]+(li-l1)/a]]), q:subst([var[1]=0,var[i]=0],q)-l1*li/a, var:append(rest(rest(var,i-1-length(var))),rest(var,i)))), display(q,var,liste)), q:sum(liste[i][1]*liste[i][2]^2,i,1,length(liste)), [liste,q])$ |