help-octave
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: sqp() with trapz()


From: Doug Stewart
Subject: Re: sqp() with trapz()
Date: Sun, 25 Aug 2019 15:11:21 -0400



On Sun, Aug 25, 2019 at 11:18 AM dg.aragones <address@hidden> wrote:
Hello all,

I have been working for a while in a code in Octave. All lines seem to work,
except for the last one, when I run an optimization using sqp(). This is the
code:

%Use Gnuplot for graphics in Octave
graphics_toolkit("gnuplot")

%Exact solution
F=@(mu) mu*log((-1-sqrt(1+mu^2))/(1-sqrt(1+mu^2)))-1; %Function for the
solver
mu=fzero(F,[0.01,1]) %Solution via Brent-Dekker-van Wijngaarden method
s=linspace(0,2,N=51); %s discretization
y=-sqrt(1+mu^2)+sqrt((s-1).^2+mu^2); %y exact values
%plot(s,y) %Exact solution

%Approximate solution
n=50; %Number of subintervals
h=2/n; %Step size
F=@(y) trapz(s,y); %Approximate objective function (Trapezoidal rule)
D=zeros(n+1,n+1);
D1=-eye(n/2)+diag(ones([n/2-1,1]),1);
D4=eye(n/2)+diag(-ones([n/2-1,1]),-1);
D(1:n/2,1:n/2)=D1;
D(end-n/2+1:end,end-n/2+1:end)=D4;
D(n/2,n/2+1)=1;
D(n/2+1,n/2+2)=1/2;
D(n/2+2,n/2+1)=-1;
D(n/2+1,n/2)=-1/2;
D=1/h*D; %Differentiation matrix
yp=@(y) transpose(D*transpose(y)); %Derivative (Finite Differences)
G=@(y) sqrt(1-yp(y).^2); %Function in the integral constraint
R=@(y) trapz(s,G(y))-1; %Approximation of the integral constraint
(Trapezoidal rule)
sqp(y,F,R) %Optimization via sequential quadratic programming

And the error it returns:

error: trapz: X and Y must have same shape
error: called from
    trapz at line 120 column 9
    sqp at line 351 column 7
    Catenary at line 28 column 1

However the trapz() function alone seems to work correctly. Does anyone know
what is happening here?

Thank you in advance for any help!



--
Sent from: https://octave.1599824.n4.nabble.com/Octave-General-f1599825.html



You have 2 different functions named F
which one do you want in the aqp function?


--
DASCertificate for 206392


reply via email to

[Prev in Thread] Current Thread [Next in Thread]