fsolve with writting automatically equations

In this code I write equations in this form in a for loop:

Dtau(i)+2*x(1)*l(i)*(s(i)^2*(cos(x(2))*cos(alpha_1(i))+sin(x(2))*sin(alpha_1(i)))

and solve the nonlinear system with fsolve for x.

When i add the definition of V_AB(i) to calculate the times for the equations it gives me this error:

Error using inlineeval

Error in inline expression ==> 56.5685424949+2*x(1)*1.0000000000*((0.0021272263)^2)*(cos(x(2))*cos(1.5707963268)+sin(x(2))*sin(1.5707963268))

0.0000000000+2*x(1)*0.7071067812*((0.0021272263)^2)*(cos(x(2))*cos(-0.7853981634)+sin(x(2))*sin(-0.7853981634)) 80.0000000000+2*x(1)*0.7071067812*((0.0021272263)^2)*(cos(x(2))*cos(0.7853981634)+sin(x(2))*sin(0.7853981634))

80.0000000000+2*x(1)*0.7071067812*((0.0021272263)^2)*(cos(x(2))*cos(0.7853981634)+sin(x(2))*sin(0.7853981634)) 0.0000000000+2*x(1)*0.7071067812*((0.0021272263)^2)*(cos(x(2))*cos(-0.7853981634)+sin(x(2))*sin(-0.7853981634))

56.5685424949+2*x(1)*1.0000000000*((0.0021272263)^2)*(cos(x(2))*cos(0.0000000000)+sin(x(2))*sin(0.0000000000))

Invalid expression. Check for missing multiplication operator, missing or unbalanced delimiters, or other syntax error. To construct matrices, use brackets instead of parentheses.

While if i calculate V_AB as a constant, without the index 'i', fsolve works.

The code i wrote is:

%defining constants
k=1.4;
R_gas=287;
T_media_vera=550; %[K]
CL_VERA=(k*R_gas*T_media_vera)^0.5;
s=1/CL_VERA;
V_vera=40; %[m/s]
Beta_vero=pi/4; %[rad]
p = rand(6,2);
q = rand(6,2);

% writing equations in a txt file:
FID = fopen('equations list.txt','w');
F=zeros(length(p),1);
for i=1:length(p)
    % lengths
    l(i)=((q(i,1)-p(i,1))^2+(q(i,2)-p(i,2))^2)^0.5;
    % angles
    alpha_1(i)=atan((q(i,2)-p(i,2))/((q(i,1)-p(i,1))));
    %velocity
    V_AB(i)=V_vera*(cos(alpha_1(i))*cos(Beta_vero)+(sin(alpha_1(i))*sin(Beta_vero)));
    % times
    tau_f(i)=l(i)/(CL_VERA+V_vera);
    tau_b(i)=l(i)/(CL_VERA-V_vera);
    %     tau_f(i)=l(i)/(CL_VERA+(V_AB(i)));  
    %     tau_b(i)=l(i)/(CL_VERA-(V_AB(i)));
    % if i add these two instead of the 2 lines above the code gives me the error!
    Dtau(i)=tau_f(i)-tau_b(i);
    %    printing
    formatSpec = '%.10f+2*x(1)*%.10f*((%.10f)^2)*(cos(x(2))*cos(%.10f)+sin(x(2))*sin(%.10f)) \n';
    F(i)=fprintf(FID,formatSpec,Dtau(i),l(i),s,alpha_1(i),alpha_1(i));
end

%importing equations and preparing them for fsolve
A=importdata('equations list.txt');
B=transpose(A);
C=cell2mat(B);

%defining problem
problem.objective = C;
problem.x0 = [0,0];
problem.solver = 'fsolve';

%setting tolerances
problem.options = optimoptions('fsolve','MaxIter', ...
    4000,'MaxFunEvals',4000,'StepTolerance',1e-16, ...
    'FunctionTolerance',1e-16,'OptimalityTolerance',1e-10);

%solving
x = fsolve(problem)

Warning: Trust-region-dogleg algorithm of FSOLVE cannot handle non-square systems; using Levenberg-Marquardt algorithm instead.

Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.

x = 1×2

58.9244 -0.5537

Does anybody know hot to solve this using

    %     tau_f(i)=l(i)/(CL_VERA+(V_AB(i)));  
    %     tau_b(i)=l(i)/(CL_VERA-(V_AB(i)));

?

Any help would be greatly appreciated.

Answer 1

[Karim]

Please try to avoid asking the same problems in many different questions. This makes it difficult for other people to track tings.

Below you can find the operational code, the easiest solution is to write a function instead of a text file. Later you can call this function to solve the problem.

%defining constants
k=1.4;
R_gas=287;
T_media_vera=550; %[K]
CL_VERA=(k*R_gas*T_media_vera)^0.5;
s=1/CL_VERA;
V_vera=40; %[m/s]
Beta_vero=pi/4; %[rad]
p = rand(6,2);
q = rand(6,2);

% writing equations in a txt file:
FID = fopen('MyFun.m','w');

% first wirte the starting line
fprintf(FID,'function F = MyFun(x)\n');

for i=1:length(p)
    % lengths
    l(i)=((q(i,1)-p(i,1))^2+(q(i,2)-p(i,2))^2)^0.5;
    % angles
    alpha_1(i)=atan((q(i,2)-p(i,2))/((q(i,1)-p(i,1))));
    %velocity
    V_AB(i)=V_vera*(cos(alpha_1(i))*cos(Beta_vero)+(sin(alpha_1(i))*sin(Beta_vero)));
    % times
    tau_f(i)=l(i)/(CL_VERA+(V_AB(i)));  
    tau_b(i)=l(i)/(CL_VERA-(V_AB(i)));
    
    Dtau(i)=tau_f(i)-tau_b(i);
    
    % printing
    formatSpec = 'F(%i) = %.10f+2*x(1)*%.10f*((%.10f)^2)*(cos(x(2))*cos(%.10f)+sin(x(2))*sin(%.10f)); \n';
    fprintf(FID,formatSpec,i,Dtau(i),l(i),s,alpha_1(i),alpha_1(i));
end

% write the 'end' section of the function
fprintf(FID,'end\n');

% close the file
fclose(FID);

% import the file as a text file just to display the result
readlines('MyFun.m')

ans = 9×1 string array
    "function F = MyFun(x)"
    "F(1) = 0.0000133197+2*x(1)*0.1821479860*((0.0021272263)^2)*(cos(x(2))*cos(-0.9887382134)+sin(x(2))*sin(-0.9887382134)); "
    "F(2) = -0.0000452168+2*x(1)*0.1749608985*((0.0021272263)^2)*(cos(x(2))*cos(0.0059349158)+sin(x(2))*sin(0.0059349158)); "
    "F(3) = -0.0001454711+2*x(1)*0.5430142149*((0.0021272263)^2)*(cos(x(2))*cos(1.5274003990)+sin(x(2))*sin(1.5274003990)); "
    "F(4) = -0.0000256197+2*x(1)*0.8365567796*((0.0021272263)^2)*(cos(x(2))*cos(-0.7007031371)+sin(x(2))*sin(-0.7007031371)); "
    "F(5) = -0.0001393017+2*x(1)*0.4375015986*((0.0021272263)^2)*(cos(x(2))*cos(1.2914281100)+sin(x(2))*sin(1.2914281100)); "
    "F(6) = -0.0001174437+2*x(1)*0.3597680897*((0.0021272263)^2)*(cos(x(2))*cos(0.3264296136)+sin(x(2))*sin(0.3264296136)); "
    "end"
    ""

%defining problem
problem.objective = @MyFun; % here just link to the function file we created
problem.x0 = [0,0];
problem.solver = 'fsolve';

%setting tolerances
problem.options = optimoptions('fsolve','MaxIter', ...
    4000,'MaxFunEvals',4000,'StepTolerance',1e-16, ...
    'FunctionTolerance',1e-16,'OptimalityTolerance',1e-10,...
    'Algorithm','levenberg-marquardt');

%solving
x = fsolve(problem);

Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.

% display the solution
x

x = 1×2

40.1940 7.0687

[anto]

OK, thanks.

But as you did, if i got this right, we solve the equations one by one.

But what if i wanted to solve them all together as system of equations?

Shouldn't we remove the last for cycle?

[Karim]

then you should write a .m file and not a .txt file

i.e. you would need to create a function that looks like the example below... i will update the answer accordingly

function F = MyFun(x)
    F(1) = -0.0003646474+2*x(1)*1.0000000000*(0.0021272263)^2*(cos(x(2))*cos(1.5707963268)+sin(x(2))*sin(1.5707963268));
    F(2) = -0.0002578447+2*x(1)*0.7071067812*(0.0021272263)^2*(cos(x(2))*cos(-0.7853981634)+sin(x(2))*sin(-0.7853981634));
    F(3) = -0.0002578447+2*x(1)*0.7071067812*(0.0021272263)^2*(cos(x(2))*cos(0.7853981634)+sin(x(2))*sin(0.7853981634));
    F(4) = -0.0002578447+2*x(1)*0.7071067812*(0.0021272263)^2*(cos(x(2))*cos(0.7853981634)+sin(x(2))*sin(0.7853981634));
    F(5) = -0.0002578447+2*x(1)*0.7071067812*(0.0021272263)^2*(cos(x(2))*cos(-0.7853981634)+sin(x(2))*sin(-0.7853981634));
    F(6) = -0.0003646474+2*x(1)*1.0000000000*(0.0021272263)^2*(cos(x(2))*cos(0.0000000000)+sin(x(2))*sin(0.0000000000));
end

[anto]

Edit: ok i saw the edit. by putting the .m file it works!!! Thank you so much you saved me

[Walter Roberson]

F = arrayfun(@str2func, "@(x)" + readlines('MyFun.m'))

Answer 2

[Walter Roberson]

    F(i) = str2func("@(x)" + sprintf(FID,formatSpec,Dtau(i),l(i),s,alpha_1(i),alpha_1(i)));

and leave out the file I/O

Notice that it is sprintf() not fprintf() here. The return value from fprintf() is the number of items converted, not a string or character vector.