Using 'for' loop in integral correctly
7 views (last 30 days)
Show older comments
Hello! I need to solve an integral which is definite and has an unknown parameter called 'B'. Its limits are t1(i) and t2(i). Because of i=1:3, I need to use the for loop. I'm trying to do it, but it's giving me results with long numbers and I felt like I'm doing something wrong. I'm sending the picture of integral I suppose to solve.
And here is the code in matlab I wrote for that integral:
for i=1:3
syms B t
y(i)=int((1-(B/t)+(1/2)*(B^2/t^2)-(1/6)*(B^3/t^3))*(a(i)+(b(i)*t)+(c(i)*t^2)),t1(i),t2(i))
end
y=y'
I tried to have the results as 'y(i)'. Because I will need yi(1) and yi(2) values for determining 'B' at the next step. I am expecting your help.
0 Comments
Accepted Answer
Alan Stevens
on 12 Aug 2020
Just define y as a function of t, B, a, b, and c, wihout regard to i. Then call the function with the appropriate parameters when needed. ie. something like:
y = @(t,B,a,b,c) int(1- B/t +(1/2)*(B/t)^2 + ... etc.);
10 Comments
Alan Stevens
on 3 Sep 2020
Ok, here's a version that allows you to continue the calculations easily:
a = [1.2235706, 1.2010072, 1.188249];
b = [-11.0411682, -9.767716, -9.0103308];
c = [23.5212358, 19.3560454, 16.8841345];
t1 = [0.0212069, 0.0214943, 0.0217816];
t2 = [0.1351592, 0.1773025, 0.2120814];
p = [1.6586, 3.6583, 4.514]*10^-5;
rho = [0.9956, 0.9922, 0.9880]*10^6;
T = [303, 313, 323];
betaT = [0.2688823, 0.271018, 0.2736945];
k = p(1)/p(2);
h1 = rho(1)*T(1)/betaT(1);
h2 = rho(2)*T(2)/betaT(2);
h = h1/h2;
B0 = 0.1; % Initial guess
B = fzero(@rfn,B0,[],a,b,c,t1,t2,k,h);
disp(B)
DI1 = yfnintegral(t2(1),B,a(1),b(1),c(1)) - yfnintegral(t1(1),B,a(1),b(1),c(1));
DI2 = yfnintegral(t2(2),B,a(2),b(2),c(2) ) - yfnintegral(t1(2),B,a(2),b(2),c(2) );
% Insert the other calculations here before the function definitions
function r = rfn(B,a,b,c,t1,t2,k,h)
DI1 = yfnintegral(t2(1),B,a(1),b(1),c(1)) - yfnintegral(t1(1),B,a(1),b(1),c(1));
DI2 = yfnintegral(t2(2),B,a(2),b(2),c(2) ) - yfnintegral(t1(2),B,a(2),b(2),c(2) );
r = h*DI1/DI2 - k;
end
function I = yfnintegral(t,B,a,b,c)
I1 = t.^2*(b - B*c)/2;
I2 = -log(t)*(a*B - b*B^2/2 + c*B^3/6);
I3 = c*t.^3/3;
I4 = -( t*(a*B^2/2 -b*B^3/6) - a*B^3/12 )./t.^2;
I5 = t*(a - b*B + c*B^2/2);
I = I1 + I2 + I3 + I4 + I5;
end
Rik
on 4 Sep 2020
This answer is the solution of the problem.
More Answers (0)
See Also
Categories
Find more on Symbolic Math Toolbox in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!