2012-05-09 3 views
1

MATLAB을 사용하여 오히려 추한 적분을 계산하려고합니다. 내가 가지고있는 문제는 아주 작은 숫자 (< 10^-300)로 매우 큰 숫자 (> 10^300)를 곱하는 부분입니다. MATLAB은 0-0.0005 범위에 있어야하지만 'inf'를 반환합니다. 이것은 내가 무엇을 가지고 큰 숫자를 작은 숫자로 곱하기

besselFunction = @(u)besseli(qb,2*sqrt(lambda*(theta + mu)).*u); 
    exponentFuncion = @(u)exp(-u.*(lambda + theta + mu)); 

QB = 5, 람다 = 12, 세타 = 10, MU = 3 그리고 내가 찾으려는 U의 모든 실제 값에 대한

besselFunction(u)*exponentFunction(u) 

입니다

. 문제는 u> 28 일 때마다 'inf'로 평가된다는 것입니다. 내가 들었고, MATLAB 함수 'vpa'를 사용하려고했지만 기능을 사용하려고 할 때 잘 작동하지 않는 것 같습니다 ...

이 시점에서 어떤 팁을 주시면 감사하겠습니다!

답변

5

대수를 사용합니다.

x = Bessel function of uy = x*exp(-u) (방정식보다 간단하지만 비슷 함). log(v*w) = log(v) + log(w) 이후

log(y) = log(x) + log(exp(-u))

이이 더 나은 수치 적으로 행동한다

log(y) = log(x) - u 

에 단순화합니다.

다른 키는 이 아니며,은 큰 숫자가되고 로그를 얻기 위해 수학 함수로 전달되는 베셀 함수를 평가합니다. Bessel 함수의 로그를 직접 반환하는 자신 만의 함수를 작성하는 것이 좋습니다. Abramowitz와 Stegun과 같은 참고 문헌을보고 시도해보십시오.

+0

실제로 가장 똑똑한 방법입니다. 감사! 내 필요에 따라 큰 숫자를 처리 할 수있는 상징적 인 엔진을 사용하는 것으로 충분했다. 이것은 내가 한 방법입니다 : besselFunction = @ (u) besseli (qb, 2 * sqrt (lambda * (theta + mu)). * sym (u)); exponentFuncion = @ (u) exp (-sym (u). * (람다 + 세타 + 뮤)); Then besselFunction (u) * exponentFunction (u) 은 기호 값을 반환합니다. 배정도 표현을 원하면 double()을 사용할 수 있습니다. – Groot

0

통합을 수행하는 경우 대신 Gauss–Laguerre quadrature을 사용해보십시오. 기본 아이디어는 형식이 exp(-x)*f(x) 인 경우 0에서 inf 로의 적분은 sum(w(X).*f(X))으로 근사 할 수 있으며 여기서 X의 값은 Laguerre 다항식의 0이며 W(X)은 특정 가중치입니다 (위키 백과 문서 참조). 매우 진보 된 심슨의 규칙과 비슷합니다. 방정식에는 이미 exp(-x) 부분이 있기 때문에 특히 적합합니다.

다항식의 근원을 찾으려면 MATLAB Central에서 LaguerrePoly이라는 함수가 있으며 여기에서 가중치를 계산하는 것은 매우 간단합니다.

관련 문제