For loop inside functions in fsolve

4 visualizaciones (últimos 30 días)
Meva
Meva el 21 de Ag. de 2016
Editada: Qu Cao el 24 de Ag. de 2016
Hello, I have 2 by 2 system of nonlinear equations.
I will use for loop for my time T inside these two functions.
I have this function for fsolve:
function F=torder1(x)
x_1=[0:0.01:1];
b=0.6;
dt = 0.01;
Nt = 1/dt+1;
%$ 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)});
*for nt=1:Nt
T(nt)=(nt-1)*dt;
F(nt)* = [double(x(1) - integral_result_1*T(nt)^2 -h0);
double(x(2) - integral_result_2*T(nt)^2 - theta0)];
end
Also I have my script in which the function torder1.m is nested:
x0 = [.1, .1];
options = optimoptions('fsolve','Display','iter');
dt=0.01;
Nt=1/dt+1;
[x,fval] = fsolve(@torder1,x0,options)
h = x(1)
theta = x(2)
I have to loop over my time T (very bottom in the torder1.m function) and in every time step,
I need to get my x(1) = h and x(2) = theta and related F in every time and save them to plot against T.
The time index is shown as nt.
Thanks for any help.

Respuestas (1)

Qu Cao
Qu Cao el 24 de Ag. de 2016
Editada: Qu Cao el 24 de Ag. de 2016
I recommend you to provide more details about the equations you want to solve.
You can also refer to the following documentation to see how to use 'fsolve':

Community Treasure Hunt

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

Start Hunting!

Translated by