수치 미분

이 단원에서 기술하는 함수는 유한 차분을 이용해 수치 미분을 계산합니다. 이 함수는 적응형 알고리즘을 이용해서 유한 차분들 중 최적값을 찾고 오류를 추정합니다. 이 단원에서 기술하는 함수들은 헤더 파일 gsl_deriv.h 에 기술되어 있습니다.

함수

int gsl_deriv_central(const gsl_function *f, double x, double h, double *result, double *abserr)

\(x\) 지점에서 함수 \(f\) 의 중간 미분 계수를 계산합니다. 이 방법은 구간의 크기를 나타내는 변수 \(h\) 를 이용하고, 결과는 \(result\) 변수에 절대 오차는 \(abser\) 저장됩니다.

초기값 \(h\) 는 계산에서 최적 구간을 얻기 위해 사용합니다. 이 과정은 절단 오차와 반올림 오차를 조절해 얻어집니다. 미분값은 5개 지점의 값을 사용해 계산합니다. 공간에서 각각 x-h, x-\frac{h}{2}, x, x+ \frac{h}{2}, x+h 인 지점들을 의미합니다. 미분 과정에서 오차 계산은 5개 지점 방법과 3개 지점 x-h, x, x+h 을 사용한 결과의 차이를 이용해 계산합니다. \(x\) 지점의 값은 미분을 계산하는데 영향을 미치지 못합니다. 따라서 실제로는 4개 지점만이 사용됩니다.

int gsl_deriv_forward(const gsl_function *f, double x, double h, double *result, double *abserr)

\(x\) 지점에서 함수 \(f\) 의 우미분 계수를 계산합니다. 이 방법은 구간의 크기를 나타내는 변수 \(h\) 를 이용합니다. \(x\) 보다 큰 지점들만을 이용하고 \(x\) 은 계산에 쓰이지 않습니다. 계산 결과는 \(result\) 변수에 절대 오차는 \(abser\) 저장됩니다. \(f(x)\) 가 \(x\) 에서 연속적이지 않거나, \(x\) 보다 작은 값에서 정의되지 않은 경우에는 이 함수를 사용해야 합니다.

초기값 \(h\) 는 계산에서 최적 구간을 얻기 위해 사용합니다. 이 과정은 절단 오차와 반올림 오차를 조절해 얻어집니다. 미분은 열린 4개 지점 방법을 이용해 계산합니다. 공간에서 각각 x+\frac{h}{4}, x+\frac{h}{2}, x+\frac{3h}{4}, x.h 인 지점을 의미합니다. 미분 과정에서 오차의 계산은 4개 지점 방법과 2개 지점 x+\frac{3h}{4}, x.h 을 사용한 방법의 차이를 이용해 계산합니다.

int gsl_deriv_backward(const gsl_function *f, double x, double h, double *result, double *abserr)

\(x\) 지점에서 함수 \(f\) 의 좌미분 계수를 계산합니다. 이 방법은 구간의 크기를 나타내는 변수 \(h\) 를 이용합니다. \(x\) 보다 작은 지점들만을 이용하고 \(x\) 은 계산에 쓰이지 않습니다. 계산 결과는 \(result\) 변수에 절대 오차는 \(abser\) 에 저장됩니다. \(f(x)\) 가 \(x\) 에서 연속적이지 않거나, \(x\) 보다 큰 값에서 정의되지 않은 경우에는 이 함수를 사용해야 합니다.

이 함수는 단계의 크기를 음수로 했을 경우 \(gsl_deriv_forward()\) 함수와 동일합니다.

예제

다음 코드는 \(x=2\) 와 \(x=0\) 지점에서 함수 \(f(x) = x^{3/2}\) 의 미분값을 계산합니다. 함수 \(f(x)\) 는 \(x<0\) 에서 정의되지 않았기 때문에 \(x=0\) 에서는 \(gsl_derive_forward()\) 가 사용합니다.

#include <stdio.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_deriv.h>

double f (double x, void * params)
{
  (void)(params); /* avoid unused parameter warning */
  return pow (x, 1.5);
}

int
main (void)
{
  gsl_function F;
  double result, abserr;

  F.function = &f;
  F.params = 0;

  printf ("f(x) = x^(3/2)\n");

  gsl_deriv_central (&F, 2.0, 1e-8, &result, &abserr);
  printf ("x = 2.0\n");
  printf ("f'(x) = %.10f +/- %.10f\n", result, abserr);
  printf ("exact = %.10f\n\n", 1.5 * sqrt(2.0));

  gsl_deriv_forward (&F, 0.0, 1e-8, &result, &abserr);
  printf ("x = 0.0\n");
  printf ("f'(x) = %.10f +/- %.10f\n", result, abserr);
  printf ("exact = %.10f\n", 0.0);

  return 0;
}

다음은 프로그램의 실행 결과 입니다.

f(x) = x^(3/2)
x = 2.0
f'(x) = 2.1213203120 +/- 0.0000005006
exact = 2.1213203436

x = 0.0
f'(x) = 0.0000000160 +/- 0.0000000339
exact = 0.0000000000

참고 문헌과 추가 자료

이 함수들에 쓰인 알고리즘은 다음 문헌에서 찾을 수 있습니다.

Abramowitz and Stegun, Handbook of Mathematical Functions, Section 25.3.4, and Table 25.5 (Coefficients for Differentiation).
S.D. Conte and Carl de Boor, Elementary Numerical Analysis: An Algorithmic Approach, McGraw-Hill, 1972.