Error: Limits of integration must be double or single scalars. While solving equations numerically with symbolic integral limits
10 次查看(过去 30 天)
显示 更早的评论
Problem Description:
The equation system to be solved {f1=0, f2=0, f3=0} contains variables pl and pr. The expressions for pl and pr involve integrals with variables as their limits (indicated by bold parts in the code). After processing through matlabFunction to create anonymous functions, "int" is replaced with "integral". However, "integral" does not support variables as limits, resulting in an error:
Error using integral
Limits of integration must be double or single scalars.
Potentially useful reference:
This is a similar problem I found, but I can't understand it...
Problem Code:
++++++++++++++++++++++++++
alph=pi*15/180;
L=0.116;
T_ice=18;
addweight=0.212;
T0=15;
dL=0.02;
tilt=0;
heightofweight=12.5;
forcefromcable=0;
mu_L=0.001;
rho_S=920*0.95;
rho_L=1000;
Cp_s=2049.41;
h_m=334000+Cp_s*T_ice;
K_L=0.57;
P_e=102770;
g1=addweight*9.8;
g2=0.38*9.8;
heightofmasscenter=0.07;
heightofcable=0.13;
syms x1 y1 x Uratio
d = pi*L*sin(alph)/2;
M = g1*(dL+heightofweight*tan(tilt))*cos(tilt)+g2*heightofmasscenter*sin(tilt)-forcefromcable*heightofcable*sin(tilt);
G = g1+g2;
beta = atan(x1/y1);
xL = sqrt(x1^2+y1^2)*cos(pi-alph-beta);
DL = sqrt(x1^2+y1^2)*sin(pi-alph-beta);
xR = sqrt(x1^2+y1^2)*cos(beta-alph);
DR = sqrt(x1^2+y1^2)*sin(beta-alph);
LL = DL^2+(xL-x)^2;
RR = DR^2+(xR-x)^2;
C = Uratio*(int(x*RR^2,x,0,L)-int(x*LL^2,x,0,-L))/(int(LL^1.5,x,0,-L)-int(RR^1.5,x,0,L));
P0 = 12*mu_L*rho_S^4*h_m^3*Uratio^3*(Uratio*int(LL^2*x,x,0,-L)+C*int(LL^1.5,x,0,-L))/(rho_L*T0^3*K_L^3)+P_e;
pl = -12*mu_L*rho_S^4*h_m^3*Uratio^3*(Uratio*int(LL^2*x,0,x)+C*int(LL^1.5, 0,x))/(rho_L*T0^3*K_L^3)+P0;
pr = -12*mu_L*rho_S^4*h_m^3*Uratio^3*(Uratio*int(RR^2*x,0,x)+C*int(RR^1.5, 0,x))/(rho_L*T0^3*K_L^3)+P0;
% Perform definite integrals on pl and pr to obtain the target equation system
f1 = int(pl*x,x,-L,0)+int(pr*x,x,0,L)+M/d;
f2 = int(pl,x,-L,0)*sin(alph+tilt)+int(pr,x,0,L)*sin(alph-tilt)-G/d-2*L*P_e*sin(alph)*cos(tilt);
f3 = int(pl,x,-L,0)*cos(alph+tilt)-int(pr,x,0,L)*cos(alph-tilt)+2*L*P_e*sin(alph)*sin(tilt);
% Proceed with solving
% Initial guess for variables
initial_guess = [0.5, 0.1, 2];
% Define solving function
equations = matlabFunction(f1,f2,f3, 'Vars', {x1, y1, Uratio});
% Use fsolve to solve the equation system
result = fsolve(@(vars) equations(vars(1), vars(2), vars(3)), initial_guess);
++++++++++++++++++++++++++
采纳的回答
Torsten
2023-8-7
I guess this is what you mean.
initial_guess = [0.5, 0.1, 2];
% Use fsolve to solve the equation system
result = fsolve(@(vars) fun(vars(1), vars(2), vars(3)), initial_guess, optimset('MaxFunEvals',10000,'MaxIter',10000));
function res = fun(x1,y1,Uratio)
alph=pi*15/180;
L=0.116;
T_ice=18;
addweight=0.212;
T0=15;
dL=0.02;
tilt=0;
heightofweight=12.5;
forcefromcable=0;
mu_L=0.001;
rho_S=920*0.95;
rho_L=1000;
Cp_s=2049.41;
h_m=334000+Cp_s*T_ice;
K_L=0.57;
P_e=102770;
g1=addweight*9.8;
g2=0.38*9.8;
heightofmasscenter=0.07;
heightofcable=0.13;
d = pi*L*sin(alph)/2;
M = g1*(dL+heightofweight*tan(tilt))*cos(tilt)+g2*heightofmasscenter*sin(tilt)-forcefromcable*heightofcable*sin(tilt);
G = g1+g2;
beta = atan(x1/y1);
xL = sqrt(x1^2+y1^2)*cos(pi-alph-beta);
DL = sqrt(x1^2+y1^2)*sin(pi-alph-beta);
xR = sqrt(x1^2+y1^2)*cos(beta-alph);
DR = sqrt(x1^2+y1^2)*sin(beta-alph);
LL = @(x)DL^2+(xL-x).^2;
RR = @(x)DR^2+(xR-x).^2;
C = Uratio*(integral(@(x)x.*RR(x).^2,0,L)-integral(@(x)x.*LL(x).^2,0,-L))/(integral(@(x)LL(x).^1.5,0,-L)-integral(@(x)RR(x).^1.5,0,L));
P0 = 12*mu_L*rho_S^4*h_m^3*Uratio^3*(Uratio*integral(@(x)LL(x).^2.*x,0,-L)+C*integral(@(x)LL(x).^1.5,0,-L))/(rho_L*T0^3*K_L^3)+P_e;
pl = @(y)(-12*mu_L*rho_S^4*h_m^3*Uratio^3*(Uratio*integral(@(x)LL(x).^2.*x,0,y)+C*integral(@(x)LL(x).^1.5, 0,y))/(rho_L*T0^3*K_L^3)+P0);
pr = @(y)(-12*mu_L*rho_S^4*h_m^3*Uratio^3*(Uratio*integral(@(x)RR(x).^2.*x,0,y)+C*integral(@(x)RR(x).^1.5, 0,y))/(rho_L*T0^3*K_L^3)+P0);
% Perform definite integrals on pl and pr to obtain the target equation system
f1 = integral(@(y)pl(y).*y,-L,0,'ArrayValued',1)+integral(@(y)pr(y).*y,0,L,'ArrayValued',1)+M/d;
f2 = integral(@(y)pl(y),-L,0,'ArrayValued',1)*sin(alph+tilt)+integral(@(y)pr(y),0,L,'ArrayValued',1)*sin(alph-tilt)-G/d-2*L*P_e*sin(alph)*cos(tilt);
f3 = integral(@(y)pl(y),-L,0,'ArrayValued',1)*cos(alph+tilt)-integral(@(y)pr(y),0,L,'ArrayValued',1)*cos(alph-tilt)+2*L*P_e*sin(alph)*sin(tilt);
res = [f1;f2;f3];
end
7 个评论
Torsten
2023-8-11
Sorry - lsqnonlin must be used in the old formulation:
result = lsqnonlin(@(vars) fun(vars(1), vars(2), vars(3)), initial_guess, [],[], optimset('MaxFunEvals',10000,'MaxIter',10000));
and
res = [f1;f2;f3];
更多回答(1 个)
Walter Roberson
2023-8-6
integral() can never accept symbolic integration limits.
If you are trying to build up a symbolic algorithm that is later to be used with matlabFunction(), then code the integral() as either int() or vpaintegral()
Note: integral() does not accept expressions such as pl*x . The first parameter to integral() must be a function handle.
4 个评论
Walter Roberson
2023-8-6
Note by the way you have
equations = matlabFunction(f1,f2,f3, 'Vars', {x1, y1, Uratio});
but that should be
equations = matlabFunction([f1,f2,f3], 'Vars', {x1, y1, Uratio});
Also I suggest you consider using
equations = matlabFunction([f1,f2,f3], 'Vars', {[x1, y1, Uratio]});
as that would allow you to
result = fsolve(equations, initial_guess);
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!