Borrar filtros
Borrar filtros

find a zero of a two-variable function

41 visualizaciones (últimos 30 días)
Douw Gerbrand
Douw Gerbrand el 20 de Abr. de 2024
Movida: Bruno Luong el 22 de Abr. de 2024
my main.m is:
clear
clc
close all
%
X0 = [9.609 , 32.288]; %initial value close to zero of function
%
f = Trajectory(X0);
The script Trajectory.m is a messy think that returns a value of f. I want to find values of tht two input variables (X0) is the initial vbalues that give f=0
  8 comentarios
Sam Chak
Sam Chak el 21 de Abr. de 2024
Is 'Trajectory' a script or a function? If it is a function, what is the output of the function?
function output = Trajectory()
output = ...;
end
Douw Gerbrand
Douw Gerbrand el 22 de Abr. de 2024
Hi Sam:
Oh wow! I have not heard the term Brachistochrone since my second year university course on applied mathematics. Yes, in some senses it is that, but with additional complications. My forces include gravity, but also two aerodunamic forces, parallel to motion and perpendicular to motion. This m,akes the script "Trajecrtory" quite complicated, and non-linear. Anyway, Bruno spotted my error and I now have a working code.
T^hanks to all!
Douw

Iniciar sesión para comentar.

Respuestas (3)

Bruno Luong
Bruno Luong el 20 de Abr. de 2024
Editada: Bruno Luong el 20 de Abr. de 2024
If you are not unlucky
Trajectory = @(x) sum(x.^2,2)-34^2-1; % Test function, use use rather own
X0 = [9.609 , 32.288]; %initial value close to zero of function
f1 = @(dx1) Trajectory(X0+[dx1 0]);
dx1 = fzero(f1,0);
X = X0+[dx1 0];
disp(X)
10.6998 32.2880
Trajectory(X)
ans = 0
fimplicit(@(x,y) Trajectory([x(:) y(:)]).')
hold on
plot(X(1),X(2),'or')
  2 comentarios
Douw Gerbrand
Douw Gerbrand el 20 de Abr. de 2024
Hi Bruno:
My "Trajectory" is a fairly complicated thing. It involves numerical solution of a pair of non-linear ODE's. I believe it would therefore not be possible to specify it using an "@" as you suggest. Is there a way to use your approach but not using the "@"?
Bruno Luong
Bruno Luong el 20 de Abr. de 2024
Editada: Bruno Luong el 20 de Abr. de 2024
As I said I gave a simple Trajectory as example, you don't have to use it but plug in your complicated code.on this
X0 = [9.609 , 32.288]; %initial value close to zero of function
f1 = @(dx1) Trajectory(X0+[dx1 0]);
dx1 = fzero(f1,0);
X = X0+[dx1 0];
Trajectory(X)
disp(X)

Iniciar sesión para comentar.


John D'Errico
John D'Errico el 20 de Abr. de 2024
A two variable function is a surface. (@Sam Chak said this, but did not take the idea to the point of completion.)
But if a surface has any solution, where it becomes zero, then almost always it will have infinitely many solutions. In fact, the tool which draws the zeros of such a function is usually called contour. That is, you want to effectively find a contour plot.
Yes, a contour plot usually is drawn with many contours at different levels, but a contour plot is what you want, all the same. You want to see only ONE contour drawn, at z == 0. Again, a contour plot can be thought of as the set of values (x,y), suze that z(x,y) is a level surface, so z(x,y) is equal to some given constant, here zero.
For an example, I'll use a symbolic function, as there are nice tools we can use. But a function handle will also work, or just an m-file. ANY function of two independent variables will work. Of course, you want it to be a continuous function, as if not, things will get very messy.
syms x y
z = x^3 - y^3 + x^2*y^2 - 3
z = 
again, we cannot simply solve for a zero of that function, as there will be infintiely many such zeros, if any exist. We might do this:
H = fcontour(z);
H.LevelList = 0;
grid on
The curve drawn in cyan is the locus of points that solve the problem z(x,y)==0, so the level set of z, at z==0. As I said, not one solution. The problem is, a tool like fsolve will produce only one solution, it will not understand there are infinitely many such solutions, just finding only one point on one of those curves, and the solution it does find will be a function of the starting value you give it.
And of course, finding the "equation" of those curves here will be a very complicated looking thing, and that is for a rather simple problem to write. We can actually find at least a set of points that will approximate those curves, using the contourc function.
So what is it you want to see? Just a plot of the solutions? Then use tools like fcontour, or contour. I could also have used fimplicit. If your goal is a set of solutions, then we could use contourc. For example...
Xv = linspace(-5,5);
Yv = linspace(-5,5);
[X,Y] = meshgrid(Xv,Yv);
% Note my careful use of the dotted operators here to vectorize the code
% and make it properly evaluated for arrays X and Y.
zfun = @(x,y) x.^3 - y.^3 + x.^2.*y.^2 - 3;
Z = zfun(X,Y);
xy = contourc(Xv,Yv,Z,[0 0])
xy = 2x212
0 -2.3787 -2.3737 -2.3603 -2.3420 -2.3238 -2.3058 -2.2881 -2.2727 -2.2706 -2.2528 -2.2354 -2.2183 -2.2017 -2.1856 -2.1717 -2.1700 -2.1547 -2.1402 -2.1267 -2.1143 -2.1032 -2.0939 -2.0865 -2.0817 -2.0800 -2.0825 -2.0906 -2.1062 -2.1328 64.0000 5.0000 4.9712 4.8990 4.7980 4.6970 4.5960 4.4949 4.4054 4.3939 4.2929 4.1919 4.0909 3.9899 3.8889 3.7978 3.7879 3.6869 3.5859 3.4848 3.3838 3.2828 3.1818 3.0808 2.9798 2.8788 2.7778 2.6768 2.5758 2.4747
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
These are the solutions, split into two curves. Unfortunately contourc puts them into one array. Looking at the first column, we see [0;64]. That tells us the first contour had z==0, and there were 64 points found on that curve. The second branch of solutions had 146 points on it that were found.
So, if you want a set of solutions, then the simple approach is to use contourc. If you want only a plot, then just use one of the contour plotting tools.
  6 comentarios
Douw Gerbrand
Douw Gerbrand el 21 de Abr. de 2024
Hi John: I tried this:
clear
clc
close all
%
X0 = [9.609 , 32.288]; %initial value close to zero of function
%
X = fminsearch(Trajectory,X
%
In the command window I see:
Undefined function or variable 'X0'.
Error in Trajectory (line 11)
velocity_init = X0(1);
Error in main (line 7)
X = fminsearch(Trajectory,X0)
what did do wrong?
Douw Gerbrand
Douw Gerbrand el 21 de Abr. de 2024
the line was:
X = fminsearch(Trajectory,X0)

Iniciar sesión para comentar.


Bruno Luong
Bruno Luong el 20 de Abr. de 2024
Movida: Bruno Luong el 22 de Abr. de 2024
@Douw Gerbrand If your objective function f is pd quadratic form of 2 variable X, then I believe you should minimize f rather than solve for f = 0.
Also if you use MATLAB ode solver, be careful that it can be an issue of optimizing parameters, see an example here : https://www.mathworks.com/matlabcentral/answers/2105041-using-fmincon-to-solve-objective-function-in-integral-form?s_tid=srchtitle
If you have only 2 parameters to estimate you could use fminsearch.
  2 comentarios
Bruno Luong
Bruno Luong el 21 de Abr. de 2024
Movida: Bruno Luong el 22 de Abr. de 2024
Try
X0 = [9.609 , 32.288]; %initial value close to zero of function
X = fminsearch(@(X)Trajectory(X),X0);
% or alternatively a simpler
X = fminsearch(@Trajectory,X0);
Or if it throw an error, please explain us at least the declaration signature of Trajectory and what does it suppose to do precisely with the input/output arguments.
Douw Gerbrand
Douw Gerbrand el 22 de Abr. de 2024
Movida: Bruno Luong el 22 de Abr. de 2024
Hi Bruno
You have spotted my error! I was not paying attention to how the variable vector X was being handled in my script "Trajectrory" I fixed that, andd not my code works.
Many thanks.

Iniciar sesión para comentar.

Productos


Versión

R2015b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by