Numerical integration in a loop

12 visualizaciones (últimos 30 días)
Muhammad Irfan
Muhammad Irfan el 6 de Feb. de 2014
Comentada: Mike Hosea el 14 de Feb. de 2014
Hi there,
I have a problem integrating numerically using the 'for' loops I want to integrate a function for different upper limits. The upper limit is also an independent variable in the function. you might understand the steps from my code below
I want to evaluate y for different t and x is the limits of the integration function
i=0
t=0:0.01:5;
a=0
for i=1:length(t);
x=(0:0.01:t);
y=@(x)-0.044*cos(omega_11.*x).*exp(-0.266.*(t-x)).*sin(13.29.*(t-x));
z(t)=quad(y,a,x)
end
plz help :)...thanks in advance

Respuestas (2)

Walter Roberson
Walter Roberson el 6 de Feb. de 2014
The function you feed in ("y" in this case) is using a vector in t that is not necessarily the same size as the "x" that is being passed in. The documentation indicates,
The function y = fun(x) should accept a vector argument x and return a vector result y, the integrand evaluated at each element of x.
but in your code if a single x were passed in, a vector of length length(t) would result, and if the length of the x that quad() chooses to pass in happens to be the same as length(t) then you are okay, but any other length(x) that quad chooses to pass in will generate an error in subtracting x-t
Also notice that you have
x=(0:0.01:t);
but you have created t as a vector. You are going to get a warning (or error) from trying to use a vector as the upper limit of a ":" operation. Did you perhaps want t(i) as the upper limit on x ?
But that x is not being used in the "y=" because x in that assignment refers to the anonymous parameter (x) for the @ function.
That leaves the vector x being passed in as the third parameter of quad, which is the "b" parameter. The "b" parameter is required to be a scalar, not a vector.
You have too many poorly defined behaviours to guess the code that would be appropriate for you.
You should be considering switching to integral2() anyhow.
  2 comentarios
Muhammad Irfan
Muhammad Irfan el 7 de Feb. de 2014
Hi Mr. Walter, Thanks for your reply
I am sorry about my code being poorly defined. I am very new to matlab. I will try to explain my problem qualitatively
I basically want to evaluate an integral equation numerically (either quad or trapz) using limits of 0 to t. This upper limit t is changing from 0 to 5 i.e. i want to evaluate from 0(lower limit)-t=0(upper limit), 0-t=0.01, 0-t=0.02 and so on until 0-t=5
The integration equation is
y(x)=-0.044*cos(omega_11.*x).*exp(-0.266.*(t-x)).*sin(13.29.*(t-x));
I hope the question is clear. I am in slight dilemma as to how to use the 'for' loops here.
Walter Roberson
Walter Roberson el 10 de Feb. de 2014
Yes, inside the "for i" loop, in your x and y lines where you refer to t, you should instead refer to t(i). In your next line where you assign to z(t) you should instead assign to z(i)
x=(0:0.01:t(i));
y=@(x)-0.044*cos(omega_11.*x).*exp(-0.266.*(t(i)-x)).*sin(13.29.*(t(i)-x));
z(i)=quad(y,a,x);

Iniciar sesión para comentar.


Mike Hosea
Mike Hosea el 10 de Feb. de 2014
Is this what you mean?
omega_11 = pi; % or whatever
t = 0:0.01:5;
Q = zeros(size(t)); % Preallocate Q for the sake of efficiency.
for k = 1:length(t)
y = @(x)-0.044*cos(omega_11.*x).*exp(-0.266.*(t(k)-x)).*sin(13.29.*(t(k)-x));
Q(k) = integral(y,0,t(k));
end
  2 comentarios
Muhammad Irfan
Muhammad Irfan el 11 de Feb. de 2014
Hi Mr Mike,
Thanks for your reply
Unfortunately, I am getting an error with the code provided by you
*??? Undefined function or method 'integral' for input arguments of type 'function_handle'.
Error in ==> trials_part_f at 67 Q(k) = integral(y,0,t(k));* *
Mike Hosea
Mike Hosea el 14 de Feb. de 2014
That's just a newer function you don't have in your version of MATLAB. Try quadgk() instead of integral().

Iniciar sesión para comentar.

Categorías

Más información sobre Numerical Integration and Differentiation en Help Center y File Exchange.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by