Is dynamic interpolation of input of function for 'quad' integration possible?
2 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
My task is to integrate a function given limited data, e.g.
max_energy = 1000;
quad(@(x) a.*x.*exp(-b.*x),0,max_energy)
where a, b are vectors, functions of energy. My problem is that it seems quad uses Simpson's rule iteratively until some tolerance is met, so the array size of a, b would need to change dynamically to match. Is it possible to interpolate between the data points I have to match quad? Must I do it in advance or use a different integration method?
0 Commenti
Risposta accettata
Andrew Newell
il 30 Gen 2012
Given that your data are irregularly spaced, you may have to use interpolation after all. If your points are (xdata,ydata), you could integrate as below:
pp = interp1(xdata,ydata,'pchip','pp');
f = @(x) ppval(pp,x);
max_energy = 1000;
quad(f,0,max_energy)
Just don't assume the integral is as accurate as the default tolerance.
3 Commenti
Andrew Newell
il 31 Gen 2012
The default method is linear. I don't know why you're getting the NaN, but if you use interp in the form yi = interp1(x,Y,xi), you have to keep recalculating the interpolating polynomials.
Più risposte (2)
Andrew Newell
il 30 Gen 2012
Any increase in accuracy based on interpolating a and b would probably be illusory. If your points are regularly spaced, you could use trapz:
y = a.*x.*exp(-b.*x);
yInt = trapz(x,y);
5 Commenti
Bård Skaflestad
il 30 Gen 2012
@Andrew
Thanks a lot -- I'd missed that part of the documentation. I obviously need to pay more attention...
Bård Skaflestad
il 30 Gen 2012
All of MATLAB's quadrature methods require an integrand that can be evaluated at vector inputs and return an equally sized vector result.
Do you mean to say that your a and b in some way depend on x? If so, you may have to implement a traditional function (i.e., .m) file that evaluates both a and b along with the resulting integrand.
For instance,
function y = integrand(x)
a = some_function(x);
b = some_other_function(x);
y = a .* x .* exp(-b .* x);
end
Does this help at all?
2 Commenti
Bård Skaflestad
il 30 Gen 2012
And then, of course, I forgot the |quad| call:
quad(@integrand, 0, max_energy)
Vedere anche
Prodotti
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!