C/C++ 전산수학 - 테일러 급수로 알아보는 삼각함수
여기서는 테일러 급수 전개를 통해 삼각함수들의 값을 구하는 방법에 대해서 알아보고, C언어 및 C++ 프로그램을 통해서 이를 구현해보도록 하겠습니다. 특정 지점에서의 함수의 미분값들을 알면, 무한한 차수의 다항식을 통해 함수의 형태를 파악할 수 있다는 것이 테일러 급수 전개의 요점입니다.
주어진 함수의 미분들을 통해 테일러 급수 전개를 방법에 대한 기본적인 내용들은 다음 포스팅에 더 자세하게 소개되어 있습니다.
사인 (sine) 함수의 미분은 코사인 (cosine) 함수에 따라 그 형태가 결정되고, 코사인 함수의 미분은 사인 함수와 연결됩니다. 이 점을 이용해서 삼각함수를 반복적으로 미분하면 어떻게 되는지를 쉽게 파악할 수 있겠죠.
먼저 사인 함수의 테일러 급수 전개에 대해서 알아봅시다.
\[ \begin{array}{rl} \sin{x} = & \displaystyle \sum_{k = 0}^{\infty} \frac{1}{k!} \left( \left. \frac{d^k}{dx^k} \sin{x} \right|_{x = n\pi} \right) (x - n\pi)^k \\ = & \displaystyle (-1)^n \sum_{k = 0}^{\infty} \frac{1}{k!} \mathcal{A}_{\textrm{sin},k} (x - n\pi)^k \\ \mathcal{A}_{\textrm{sin},k} = & \displaystyle \left\{ \begin{array}{ll} 1 & \textrm{if} \quad k \mod 4 = 1 \\ -1 & \textrm{if} \quad k \mod 4 = 3 \\ 0 & \textrm{otherwise} \end{array} \right. \end{array} \]
사인 함수와 코사인 함수의 테일러 급수 전개는 무한대의 수렴반경을 가지고 있기 때문에 기준점을 어디에 잡아도 테일러 급수는 수렴하게 됩니다. 다만 더 빠른 계산을 위해서 원주율의 정수배가 되는 지점을 기준점으로 잡아볼 수 있습니다. 주어진 \(x\)값에 대해 \(x = n\pi + x_{\scriptsize \textrm{mod}}\)를 만족하면서 \(x_{\scriptsize \textrm{mod}}\)의 절대값이 최소가 되는 \(n\)의 값을 찾을 수 있겠죠. 원주율의 정수배가 되는 지점에서 \( \sin{(n\pi)} = 0\) 및 \( \cos{(n\pi)} = (-1)^n \) 가 된다는 점을 이용하면 삼각함수의 미분들을 쉽게 계산할 수 있겠습니다.
참고로 \( a \mod q \) 는 modulo라는 함수로서, \(a\)를 \(q\)로 나눴을 때의 나머지를 말합니다. 예를 들면 \(13 \mod 4 = 1\) 이 되죠.
코사인 함수 역시 다음과 같이 테일러 전개가 가능합니다.
\[ \begin{array}{rl} \cos{x} = & \displaystyle \sum_{k = 0}^{\infty} \frac{1}{k!} \left( \left. \frac{d^k}{dx^k} \cos{x} \right|_{x = n\pi} \right) (x - n\pi)^k \\ = & \displaystyle (-1)^n \sum_{k = 0}^{\infty} \frac{1}{k!} \mathcal{A}_{\textrm{cos},k} (x - n\pi)^k \\ \mathcal{A}_{\textrm{cos},k} = & \displaystyle \left\{ \begin{array}{ll} 1 & \textrm{if} \quad k \mod 4 = 0 \\ -1 & \textrm{if} \quad k \mod 4 = 2 \\ 0 & \textrm{otherwise} \end{array} \right. \end{array} \]
사인 함수의 경우와 비슷하게 원주율의 정수배가 되는 지점을 기준점으로 잡고, 미분값들을 계산하여 테일러 급수 전개를 할 수 있습니다.
테일러 급수 전개를 통해 계산한 삼각 함수의 값이 C언어 수학 함수 라이브러리로부터 구한것과 일치한다는 것을 볼 수 있습니다.
테일러 급수로부터 삼각함수들을 계산하기 위해 사용한 C/C++ 프로그램을 소개합니다. 범용성을 높이기 위해, 서로 다른 계산을 수행하는 함수들을 서로 다른 소스 파일과 헤더 파일에 포함시켰습니다.
가장 먼저 할 일은 주어진 \(x\)값에 대해 \(x = n\pi + x_{\scriptsize \textrm{mod}}\)를 만족시키는 \(n\)과 \(x_{\scriptsize \textrm{mod}}\)의 값들을 구해는 것인데요. 이를 계산하는 함수인 get_modulo의 프로토타입을 선언하기 위한 헤더파일이 필요합니다.
math_common.h [다운로드]
#ifndef MATH_COMMON_H #define MATH_COMMON_H #undef __WRAP_CXX_INI #undef __WRAP_CXX_FIN #ifdef __cplusplus #define __WRAP_CXX_INI extern "C" { #define __WRAP_CXX_FIN } #else #define __WRAP_CXX_INI /* empty */ #define __WRAP_CXX_FIN /* empty */ #endif __WRAP_CXX_INI extern FILE *ptr_output_math_; extern double eps_precision_; double get_modulo(double x, double modulus, int *ptr_n); __WRAP_CXX_FIN #endif
__cplusplus라는 식별자를 통해 C언어와 C++ 소스파일에서 동시에 사용가능한 범용 헤더파일을 구성할 수 있습니다. get_modulo 함수를 C언어 소스파일에서 정의했더라도 extern "C" 코드블록 덕분에 name mangling이 일어나지 않죠. 결과적으로 C++ 소스코드에서 헤더파일을 포함시키고, 함수를 호출하는데 문제가 없습니다.
전역변수 eps_precision_ 이 선언되어 있는것을 볼 수 있는데, 이 실수형 변수는 급수 전개의 수렴여부를 결정하는데 사용될 것입니다. 삼각함수의 테일러 급수전개는 원칙적으로는 무한급수입니다만, 컴퓨터가 영원히 계산을 할 수는 없기 때문에 이를 유한급수로 가정하는 근사법이 필요한데요. 새로 추가되는 항 (term)의 크기가 기존 유한급수에 비해 무시할 정도로 작다면, 수렴한다고 가정하고 급수 전개를 중단하는 것입니다.
다음 단계는 별도의 C언어 소스파일에 get_modulo 함수를 정의하는 것입니다.
math_common.c [다운로드]
#include<stdio.h> #include<stdlib.h> #include<math.h> FILE *ptr_output_math_ = NULL; double eps_precision_ = 1e-8; double get_modulo(double x, double modulus, int *ptr_n) { int n = 0; double x_mod = x; int sign = 1; if (x_mod < 0.) { sign = -1; x_mod = -x_mod; } while (x_mod > modulus * (1. + (double)n)) { n += 1; x_mod -= modulus; } if (fabs(x_mod - modulus * (double)n) > fabs(modulus * (1. + (double)n) - x_mod)) { n += 1; x_mod -= modulus; } if (sign == -1) { n = -n; x_mod = -x_mod; } *ptr_n = n; return x_mod; }
get_modulo 함수는 두 개의 실수형 변수와 한 개의 포인터를 매개변수로 받아서 실수형 변수의 값을 리턴하는데요. \( x = n q + x_{\scriptsize \textrm{mod}} \) 에서 \( x \) 와 \( q \)가 각각 첫번째와 두번째 매개변수로 들어가고, \( n \)의 값이 저장되는 정수형 변수를 가리키는 포인터가 세번째 매개변수로 들어갑니다. 마지막으로 \( x_{\scriptsize \textrm{mod}} \)의 값을 리턴합니다. 앞에서 언급한 테일러 급수 전개에 get_modulo 함수를 사용하기 위해서 \( n q \)가 \( x \)에 가장 근접해지는 경우의 \( n \)을 선택하도록 했습니다.
급수의 수렴여부를 결정하는 변수 eps_precision_의 값은 원하는 정밀도에 따라 달라질 것입니다. 더 정밀한 계산을 원할수록 eps_precision_의 값을 작게 해야 합니다.
삼각함수들을 계산하기 위한 get_sin (사인), get_cos (코사인), get_tan (탄젠트) 함수들의 프로토타입을 별도의 헤더파일에 선언합니다.
math_taylor_tri.h [다운로드]
#ifndef MATH_TAYLOR_TRI_H #define MATH_TAYLOR_TRI_H #include"math_common.h" #undef __WRAP_CXX_INI #undef __WRAP_CXX_FIN #ifdef __cplusplus #define __WRAP_CXX_INI extern "C" { #define __WRAP_CXX_FIN } #else #define __WRAP_CXX_INI /* empty */ #define __WRAP_CXX_FIN /* empty */ #endif __WRAP_CXX_INI double get_sin(double x); double get_cos(double x); double get_tan(double x); __WRAP_CXX_FIN #endif
math_common.h 헤더파일과 마찬가지로 __cplusplus 식별자 및 extern "C" 코드블록을 통해 C언어와 C++ 소스파일 모두에서 사용할 수 있습니다.
마지막으로 남은 일은 get_sin, get_cos 함수들을 앞에서 설명한 테일러 급수 전개를 통해 정의하는 것입니다. 탄젠트 값을 리턴하는 get_tan 함수는 \( \tan{x} = \sin{x} / \cos{x} \) 라는 공식을 사용하면 되겠죠.
math_taylor_tri.c [다운로드]
#include<stdio.h> #include<stdlib.h> #include<math.h> #include"math_common.h" double get_sin(double x) { int n; double x_mod = get_modulo(x, M_PI, &n); double sign_n; if (n % 2 == 0) { sign_n = 1.; } else { sign_n = -1.; } double val_sin = 0.; int converging = 0; int k = -1; double term_add = 1.; while (converging == 0) { k += 1; if (k != 0) { term_add *= (x_mod / (double)k); } double sign_k = 0.; if (k % 4 == 1) { sign_k = 1.; } else if (k % 4 == 3) { sign_k = -1.; } else { continue; } double val_sin_prev = val_sin; val_sin += sign_k * term_add; double val_sin_now = 0.5 * (fabs(val_sin + val_sin_prev) + 1.); if (fabs(val_sin - val_sin_prev) < eps_precision_ * val_sin_now) { break; } } val_sin = sign_n * val_sin; return val_sin; } double get_cos(double x) { int n; double x_mod = get_modulo(x, M_PI, &n); double sign_n; if (n % 2 == 0) { sign_n = 1.; } else { sign_n = -1.; } double val_cos = 0.; int converging = 0; int k = -1; double term_add = 1.; while (converging == 0) { k += 1; if (k != 0) { term_add *= (x_mod / (double)k); } double sign_k = 0.; if (k % 4 == 0) { sign_k = 1.; } else if (k % 4 == 2) { sign_k = -1.; } else { continue; } double val_cos_prev = val_cos; val_cos += sign_k * term_add; double val_cos_now = 0.5 * (fabs(val_cos + val_cos_prev) + 1.); if (fabs(val_cos - val_cos_prev) < eps_precision_ * val_cos_now) { break; } } val_cos = sign_n * val_cos; return val_cos; } double get_tan(double x) { return get_sin(x) / get_cos(x); }
여러개의 \( x \) 값에 대해 \( \sin{x} \) 및 \( \cos{x} \)의 값들을 계산하고, 이를 텍스트 파일에 출력하는 main 프로그램을 작성해 볼 수 있습니다.
test_taylor_trigonometry.cpp [다운로드]
#include<stdio.h> #include<stdlib.h> /* _USE_MATH_DEFINES must be * defined for Visual Studio */ //#define _USE_MATH_DEFINES #include<math.h> #include"math_taylor_tri.h" //#define _WRITE_OUTPUT_FILE_ int main(int argc, char *argv[]) { double xmin = -2. * M_PI; double xmax = 2. * M_PI; int nbin_x = 40; double delta_x = fabs(xmax - xmin) / (double)nbin_x; #ifdef _WRITE_OUTPUT_FILE_ ptr_output_math_ = fopen("tab_taylor_trigonometry.dat", "w"); #else ptr_output_math_ = stdout; #endif if (ptr_output_math_ != NULL) { fprintf(ptr_output_math_, "# [x]"); fprintf(ptr_output_math_, " [sin(x) from taylor exp]"); fprintf(ptr_output_math_, " [sin(x) from math.h lib]"); fprintf(ptr_output_math_, " [cos(x) from taylor exp]"); fprintf(ptr_output_math_, " [cos(x) from math.h lib]"); fprintf(ptr_output_math_, "\n"); } for (int ix = 0; ix <= nbin_x; ix++) { double x_now = xmin + delta_x * (double)ix; double series_sin = get_sin(x_now); double series_cos = get_cos(x_now); if (ptr_output_math_ == NULL) { continue; } fprintf(ptr_output_math_, " %e %e %e %e %e\n", x_now, series_sin, sin(x_now), series_cos, cos(x_now)); } #ifdef _WRITE_OUTPUT_FILE_ fclose(ptr_output_math_); #endif }
이번 포스팅에서 다룬 테일러 급수 전개를 진행하는데 있어서 매우 중요한 요소 중 하나가 바로 원주율 \( \pi \)입니다. 여기서는 수학 함수 라이브러리에 있는 값을 사용했습니다만, 정적분을 통해 그 값을 계산해볼 수 있습니다.
댓글
댓글 쓰기