# Matlab limitation in fsolve using function input

14 views (last 30 days)
Meva on 22 Aug 2016
Commented: Meva on 31 Aug 2016
Hello,
I tried to loop for time value (T) inside my fsolve, but fsolve is pretty unforgiving.
The time loop does not seem working.
When I plot, it gives the same values (h=x(1) and theta=x(2) does not change over time which should change)!
Please see the the script that uses for loop for time (T). T is input for fsolve. :
x0 = [.1, .1];
options = optimoptions('fsolve','Display','iter');
dt=0.01;
Nt=1/dt+1;
Tarray = [0:dt:1];
T = 0;
for nt=1:Nt
[x,fval] = fsolve(@torder1,x0,options,T)
T=T+dt;
h(nt)=x(1);
theta(nt) = x(2);
plot(Tarray,h,'*')
hold on
plot(Tarray,theta,'+')
end
and the function for fsolve:
function F=torder1(x,T)
x_1=[0:0.01:1];
b=0.6;
%\$ sol(1) = h; sol(2) =theta;
clear x_1;
syms x_1 h theta kappa
f_1(x_1,h,theta,kappa) = 1/2*(1-( (h+(1-b)*theta)^2/(h+(x_1-b)*theta-kappa*x_1*(1-x_1))^2 ));
f_2(x_1,h,theta,kappa) = (x_1-b)/2*( 1-( (h+(1-b)*theta)^2/(h+(x_1-b)*theta-kappa*x_1*(1-x_1))^2 ));
kappa =1;
f_11 = 1-( (h+(x_1-b)*theta)^2/(h+(x_1-b)*theta-1*x_1*(1-x_1))^2 );
f_21 = (x_1-b)/2*( 1-( (h+(1-b)*theta)^2/(h+(x_1-b)*theta-x_1*(1-x_1))^2 ));
fint_1 = int(f_11, x_1);
fint_2 = int(f_21, x_1);
x_1=1;
upper_1=subs(fint_1);
upper_2=subs(fint_2);
clear x_1;
x_1=0;
lower_1=subs(fint_1);
lower_2=subs(fint_2);
clear x_1;
integral_result_1old=upper_1-lower_1;
integral_result_2old=upper_2-lower_2;
h0 = kappa *b*(1-b);
theta0 = kappa*(1-2*b);
integral_result_1 = subs(integral_result_1old, {h, theta}, {x(1), x(2)});
integral_result_2 = subs(integral_result_2old, {h, theta}, {x(1), x(2)});
F = [double(x(1) - integral_result_1*T^2 -h0);
double(x(2) - integral_result_2*T^2 - theta0)];

Walter Roberson on 22 Aug 2016
Edited: Walter Roberson on 22 Aug 2016
fsolve() is for real values, but solutions to your equations are strictly complex, except at T = 0.
You might be getting false results from fsolve(), with it either giving up or finding something that appears to come out within constraints.
For example, if you
fsolve(@(x) x^2+1, rand)
then you will get a small negative real-valued answer that MATLAB finds to be within the tolerances, when the right answer should be something close to sqrt(-1)

Meva on 25 Aug 2016
If I do omit double expression inside F function, and then type :
X = sym('x', [1,2]);
T = 1/100;
F = torder1(X,T);
sols = vpasolve(F, X);
sols.x1
sols.x2
This gives the following error :
CAT arguments dimensions are not consistent.
Alternatively, If I omit double expression inside F function, and then type :
X = sym('x', [1,2]);
T = 1/100;
F = torder1(X,T);
sols = vpa(F);
sols.x1
sols.x2
It still gives the CAT error as above.
So I conclude that there is a conflict between prescribing your function arguments in the type of double, and VPA and/or vpasolve.
Walter Roberson on 27 Aug 2016
I sleep. I have drive failures. I have network failures. I have appointments.
Meva on 31 Aug 2016
A sincere apologise Walter. Sorry.

Alan Weiss on 22 Aug 2016
I think that you need to replace your line
[x,fval] = fsolve(@torder1,x0,options,T)
with
[x,fval] = fsolve(@(x)torder1(x,T),x0,options)
Alan Weiss
MATLAB mathematical toolbox documentation

Meva on 22 Aug 2016
Nope!
[x,fval] = fsolve(@torder1(x,T),x0,options)
is invalid MATLAB syntex.
John D'Errico on 22 Aug 2016
But that is not what Alan suggested. There is a difference between these lines:
[x,fval] = fsolve(@torder1(x,T),x0,options)
[x,fval] = fsolve(@(x)torder1(x,T),x0,options)
Alan suggested the second, but you then tried the first.
Meva on 22 Aug 2016
You are right. However, it does give me again incorrect plot as attached. It does not understand my loop!