Hi,
I'm trying to implement solution of a nonlinear system, in which i'd like to use a scatteredInterpolant to calculate some values. 9 equations.
I have 3 vectors containing the data in which to construct the interpolant: Q, H and P_idr_gr1 , that is power of turbine vs. discharge and head.
The interpolant used from the command window works just fine.
When i try to implement in the nonlinear system i get the errors:
Index in position 1 is invalid. Array indices must be positive integers or logical values. (or: Index in position 1 exceeds array bounds, depending on x0)
Error in nlsystem>mysystem (line 55)
S(7) = x(8)*x(4)-pe1*1000/gamma/etac(x(8),x(4));
Error in fsolve (line 242)
fuser = feval(funfcn{3},x,varargin{:});
Error in nlsystem (line 41)
[x1,fval1,exitflag1] = fsolve(@mysystem,x0,options);
Caused by:
Failure in initial objective function evaluation. FSOLVE cannot continue.
I'd really appreciate to understand what i'm, doing wrong !
Thanks,
Cristiano
% my nonlinear system using scatteredInterpolant
clear allclc;global Hm pe1 pe2 pe3 etae2 etae3 qmax1 qmax2 qmax3%% dati iniziali
load '190607 collinare USBR fig 35 rendimento turbina.mat'Hm = 485.90;pe1 = 10.0;pe2 = 1.40;pe3 = 1.40;etae_2 = [-0.0131 0.042 0.0138 0.7827]; etae2 = polyval(etae_2, pe2);etae_3 = [-0.0131 0.042 0.0138 0.7827]; etae3 = polyval(etae_3, pe3);% here i create surface of hill chart
Q100 = 27.9;H100 = 81;Q = Q100.*cUSBRfig35.rel_q;H = H100.*cUSBRfig35.rel_h;P_idr_gr1 = 9.8045*Q.*H;% remove generator losses and calculate efficiency eta
P_erog_gr1 = P_idr_gr1.*cUSBRfig35.eta_t - (0.0089.*P_idr_gr1+291.77);eta_erog_gr1 = P_erog_gr1./P_idr_gr1;% calculate interpolated surface
etac = scatteredInterpolant(Q,H,eta_erog_gr1,'linear');% some initial condition for fsolve
x0 = [ Hm Hm Hm (Hm-387.29) (Hm-387.29) 387.29 0.0 0.0 0.0 ];%% lancio solutore
options = optimoptions('fsolve','Display','iter', 'FunctionTolerance', 1e-6);[x1,fval1,exitflag1] = fsolve(@mysystem,x0,options);function S = mysystem(x) global Hm pe1 pe2 pe3 etae2 etae3 etac gamma = 9.8045; S(1) = Hm -x(1) -0.018588*x(7)^2; S(2) = x(1)-x(2)-0.000541*x(8)^2; S(3) = x(1)-x(3)-0.007589*x(9)^2; S(4) = x(4)-x(2)+x(6); S(5) = x(5)-x(3)+x(6); S(6) = x(6) - (0.0003*x(7)^3 - 0.0183*x(7)^2 + 0.4374*x(7) + 387.29); S(7) = x(8)*x(4)-pe1*1000/gamma/etac(x(8),x(4)); S(8) = x(9)*x(5)-(pe2*1000/gamma/etae2 + pe3*1000/gamma/etae3); S(9) = x(7)-x(8)-x(9);end
Best Answer