Fsolve for17 simultaneous equations
13 views (last 30 days)
Show older comments
Hello,
I am attempting to set up an optimization code for solving a fluid dynamics problem. Since there are many unique variables, I wanted to know if there was a way to define each variable as part of a vector or array? I would like to be able to write my equations in terms of their native variables instead of x(0),x(1),x(2)...x(n).
Example: A = p(1), B = p(2), C = p(3)...., Q = p(17)
I realize that I there is likely a much better way, but this is the first problem of this magnitude that I have attempted. I have only used fsolve for single variable iterations prior to this and was unable to find an answer on this forum or Google. Likely due to the fact that I do not know the correct question to ask.
Thanks so much for the help.
%function [F] = Fluids_Project_2a
format long;
options = optimset('TolFun', 1e-4, 'TolX', 1e-4, 'MAxIter',1e5,'MaxFunEvals',1e5);
%Variables
p{1}=Re1;
p{2}=Re2;
p{3}=Re3;
p{4}=u1;
p{5}=u2;
p{6}=u3;
p{7}=fF1;
p{8}=fF2;
p{9}=fF3;
p{10}=F1;
p{11}=F2;
p{12}=F3;
p{13}=P1;
p{14}=P2;
p{15}=Q1;
p{16}=Q2;
p{17}=Q3;
%%Equations
%Effective Length/Diameter
LD1 = L1/D1 + E_L1*Eb90 + Ga_L1*Gate+ Contraction2to1;
LD2 = L2/D2 + E_L2*Eb90 + Ga_L2*Gate + Gl_L2*Globe + Contraction2to1;
LD3 = L3/D3 + E_L3*Eb45 + Ga_L3*Gate + Gl_L3*Globe;
%
F = {@(p) (rho*p(4)*D1/mu - p(1));
(Q1*4/pi/D1^2 - u1);
((-1.737*log(0.269*epsilon/D1 - 2.185/p(1)*log(0.269*epsilon/D1 + 14.5/p(1))))^-2 - fF1);
(2*fF1*u1^2*(LD1) - F1);
%
%Pipe 2
(rho*u2*D2/mu - p(2)); %Re2
(Q2*4/pi/D2^2 - u2); %u2
((-1.737*log(0.269*epsilon/D2 - 2.185/p(2)*log(0.269*epsilon/D2 + 14.5/p(2))))^-2 - fF2); %Re2
(2*fF2*u2^2*(LD2) - F2); %F2
%
%Pipe 3
(rho*u3*D3/mu - p(3)); %Re
(Q3*4/pi/D3^2 -u3); %u3
((-1.737*log(0.269*epsilon/D3 - 2.185/p(3)*log(0.269*epsilon/D3 + 14.5/p(3))))^-2 - fF3); %fF3
(2*fF3*u3^2*(LD3) - F3); %F3
%
%Energy Balance
((Z_L2 - Z_L1) + AA*(P2 - P1 - 10)/(rho*g) + BB*F1/g); %Node 1 -> 2
(AA*(P3 - P2)/(rho*g) + BB*F2/g); %Node 2 -> 3
((Z_L3 - Z_L2) + AA*(P4 - P2)/(rho*g) + BB*fF3*Q3^2*(LD3)/g/D3^4); %Node 2 -> 4
(P1 - (a-b*Q1^2));
%
%Mass flow
(Q1 - Q2 - Q3)};
%
p0 = [10000 10000 10000 10 10 10 .001 .001 .001 200 200 200 20 20 50 50 50];
h = fsolve(f,p0,options);
0 Comments
Accepted Answer
More Answers (0)
See Also
Categories
Find more on Function Creation in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!