수치해석 161009 Ch21. Numerical Differentiation
21.1 소개및배경 (1/2) 미분 도함수 : 독립변수에대한종속변수의변화율 y = x f ( xi + x) f ( xi ) x dy dx f ( xi + x) f ( xi ) = lim = y = f ( xi ) x 0 x 차분근사 도함수 1 차도함수 : 곡선의한점에서접선의구배
21.1 소개및배경 (2/2) 2 차도함수 : 구배가얼마나빨리변하는지의정도, 곡율 2 d y d dy 2 dx dx dx = 편도함수는한개이상의변수에의존하는함수에대해사용한다. 한개의변수를고정시키고, 한점에서함수의도함수를취한다. f f ( x + x, y) f ( x, y) = lim x x 0 x f f ( x, y + y) f ( x, y) = lim y y 0 y
21.2 고정확도미분공식 (1/4) 정확도의수준은사용하는 Taylor series 의항의개수에비례. 예를들면, 다시쓰면, 전향차분공식 (2 차이상의도함수항을무시 ) f ( xi+ 1 ) f ( xi ) f ( xi ) = + O( h) h f ( xi ) 2 f ( xi+ 1) = f ( xi ) + f ( xi ) h + h +L 2! f ( x + ) f ( x ) f ( x ) f x = h + O h h 2! i 1 i i 2 ( i ) ( ) 더정확한공식은? 2 차도함수에대한전향차분근사 f ( xi+ 2) 2 f ( xi + 1) + f ( xi ) f ( xi ) = + O( h) 2 h 위식에대입하면, f ( x ) + 4 f ( x ) 3 f ( x ) f x = + O h 2h i+ 2 i+ 1 i 2 ( i ) ( ) 2 차도함수를고려함으로써정확도를 O(h 2 ) 으로향상시켰다.
21.2 고정확도미분공식 (2/4) Forward Finite-Difference
21.2 고정확도미분공식 (3/4) Backward Finite-Difference
21.2 고정확도미분공식 (4/4) Centered Finite-Difference
예제 21.1 (1/2) Q. 예제 4.4 에서다음함수의도함수를유한차분과간격크기 h=0.25 를이용하여 x=0.5 에서구하였다. f x x x x x 4 3 2 ( ) = 0.1 0.15 0.5 0.25 + 1.2 그결과는다음표와같다. 후향 O(h) 중심 O(h 2 ) 전향 O(h) 추정값 -0.714-0.934-1.155 ε t 21.7% -2.4% -26.5% 여기서오차는정해 f'(0.5)= -0.9125 에기초한다. 고정확도공식을사용하여다시계산하라.
예제 21.1 (2/2) Sol) 필요한데이터 x x x x x = 0 f ( x ) = 1.2 i 2 i 2 = 0.25 f ( x ) = 1.1035156 i 1 i 1 i = 0.5 f ( x ) = 0.925 = 0.75 f ( x ) = 0.6363281 i+ 1 i+ 1 = 1 f ( x ) = 0.2 i+ 2 i+ 2 i O(h 2 ) 의전향차분 0.2 + 4(0.6363281) 3(0.925) f (0.5) = = 0.859375 ε t =5.82 % 2(0.25) O(h 2 ) 의후향차분 3(0.925) 4(1.1035156) + 1.2 f (0.5) = = 0.878125 ε t =3.77 % 2(0.25) O(h 4 ) 의중심차분 0.2 + 8(0.6363281) 8(1.1035156) + 1.2 f (0.5) = = 0.9125 ε t =0 % 12(0.25)
21.3 Richardson 외삽법 (1/2) 유한차분을사용할때도함수추정값을향상시키는두가지방법 -간격크기 (h) 의축소 -더많은점을포함하는고차공식의사용 Richardson 외삽법 두개의저차의도함수추정값을이용하여고차정확도의추정값을계산할수있게한다. 18장의 Richardson 외삽법은다음과같은공식을이용하여향상된적분추정값을제공 1 I = I( h2 ) + [ I( h2 ) I( h1 )] 2 ( h / h ) 1 1 2 편의상 h 2 =h 1 /2 을이용하면, 4 1 I = I( h2 ) I( h1 ) 3 3
21.3 Richardson 외삽법 (2/2) O(h 2 ) 의중심차분근사에이공식을사용하면, O(h 4 ) 의도함수추정값을구할수있다. 4 1 D = D( h2 ) D( h1 ) 3 3 O(h 4 ) 의중심차분근사에이공식을사용하면, O(h 6 ) 의도함수추정값을구할수있다. D = 16 15 D(h 2 ) 1 15 D(h 1 ) O(h 6 ) 의중심차분근사에이공식을사용하면, O(h 8 ) 의도함수추정값을구할수있다. D = 64 63 D(h 2) 1 63 D(h 1)
예제 21.2 Q. 예제 21.1 의함수와간격 h 1 =0.5, h 2 =0.25 를사용하여 x=0.5 에서의 1 차도함수를구하라. 그리고 Richardson 외삽법을사용하여향상된추정값을구하라. 정해는 -0.9125 이다. Sol) f x x x x x 4 3 2 ( ) = 0.1 0.15 0.5 0.25 + 1.2 중심차분으로 1차도함수를계산하면 0.2 1.2 D(0.5) = = 1.0 ε t = 9.6% 1 0.6363281 1.103516 D(0.25) = = 0.934375 ε t = 2.4% 0.5 Richardson 외삽법을이용한개선된추정값 4 1 D = ( 0.934375) ( 1) = 0.9125 3 3 고려하는함수가 4 차다항식이므로정확한값을도출한다.
21.4 부등간격데이터에대한도함수 부등간격데이터의도함수를계산하는한가지방법은다항식보간을수행한후도함수를구하는것이다. 예로서, 세개의점을지나는 2 차의 Lagrange 다항식을구하고, 이다항식의도함수를구하면다음과같다. f ( x)= f x 0 2x x 1 x 2 2x x 0 x 2 2x x 0 x 1 ( ) ( x 0 x 1 )( x 0 x 2 ) + f ( x 1) ( x 1 x 0 )( x 1 x 2 ) + f ( x 2) ( x 2 x 0 )( x 2 x 1 ) 장점 세점으로주어진구간내의어떤점에서도도함수를구할수있다. 주어진점들이등간격으로분포되지않아도된다. 도함수값이중심차분값의정확도와같다. 특히등간격일때 x=x 1 에서, f ( x2) f ( x0 ) 2 f ( x1 ) = + O( h ) 식 (4.25) 2h
예제 21.3 (1/2) Q. 땅속의온도구배에대한측정값은그림과같다. 흙과공기의경계에서의열플럭스는 Fourier 법칙에따라다음과같이계산된다. dt q( z = 0) = k dz z = 0 q(x)= 열플럭스 (W/m2), k= 흙의열전도계수 (=0.5 W/(m K), T= 온도 (K), z= 경계면으로부터의땅속깊이 수치미분을사용하여흙과공기의경계에서의온도구배를구하고, 이를이용하여땅속으로의열플럭스를계산하라.
예제 21.3 (2/2) Sol) 흙과공기의경계에서도함수를구하면, 2(0) 0.0125 0.0375 2(0) 0 0.0375 f 0 = 13.5 + 12 ( ) ( 0 0.0125)( 0 0.0375) ( 0.0125 0)( 0.0125 0.0375) 2(0) 0 0.0125 + 10 = 1440 + 1440 133.333 = 133.333 K / m ( 0.0375 0)( 0.0375 0.0125) 이를이용하여열플럭스를구하면, W W W q( z = 0) = 0.5 133.333 = 66.667 m K m m 2
21.5 오차를가지는데이터에대한도함수와적분 실험데이터의미분에관련된문제점은데이터에포함된오차를증폭시키는것이며, 반면에수치적분은데이터의오차를완화시킨다. 오차를가지는데이터에대해서도함수를결정하는방법은최소제곱회귀분석을사용하여데이터에매끄럽고미분가능한함수를접합한후, 그함수의도함수를구하는것이다. 데이터에포함된작은오차가 수치미분을 통해 얼마나 증폭되는지를나타내는그림 : (a) 오차가없는데이터, (b) 곡선 (a) 의수치미분결과, (c) 조금수정된데이터, 그리고 (d) 증가된변동을분명히보여주는 곡선 (c) 의 수치미분 결과. 대조적으로 반대 연산인 적분 [ 곡선 (d) 아래의면적을 취함으로써 (d) 에서 (c) 로 이동함 ] 은 오차를 줄이거나 완화시킨다.
21.6 편도함수 일차원편도함수는도함수와같은방법으로계산. f f ( x + x, y) f ( x x, y) = x 2 x f f ( x, y + y) f ( x, y y) = y 2 y 고차도함수를구하기위해서두개이상의변수에대한함수를미분할수도있다. 그결과를혼합편도함수라칭한다. f f 2 ( x + x, y) ( x x, y) f f y y = = x y x y 2 x 2 f f ( x + x, y + y) f ( x + x, y y) f ( x x, y + y) + f ( x x, y y) = x y 4 x y
21.7 MATLAB 을이용한수치미분 yi = diff(xi) xi = n 개의요소를가지는 1 차원벡터 yi = n-1 개의요소를가지는 1 차원벡터 xi 벡터내인접한요소두개값의차이 diff(y)./diff(x): 1 차도함수의유한차분근사값을계산하는데사용. fx = gradient(f) f = 원소의수가 n 인 1 차원벡터 fx = f 에기초하여계산된 n 개의차분을포함하는벡터 첫번째점은전향차분, 마지막점은후향차분, 그사이점들은중심차분을이용하여계산 편미방도계산가능
예제 21.4 (diff 을사용한미분 ) (1/2) Q. MATLAB 함수 diff 를사용하여 x=0 에서 0.8 까지다음함수를미 분하라. f(x)=0.2+25x-200x 2 +675x 3-900x 4 +400x 5 정해는 f'(x)=25-400x+2025x 2-3600x 3 +2000x 4 이다. Sol). >> f = @(x) 0.2+25*x-200*x.^2+675*x.^3-900*x.^4+400*x.^5; >> x = 0 : 0.1 : 0.8 ; >> y = f(x) ; >> diff(x) 벡터내인접한요소두개값의차이를구한다. ans = 0.1000 0.1000 0.1000 0.1000 0.1000 0.1000 0.1000 0.1000 >>d=diff(y)./diff(x) % d=diff(f(x))/0.1; 도함수의제차분근사값을구한다. d = 10.8900-0.0100 3.1900 8.4900 8.6900 1.3900-11.0100-21.3100 >> n=length(x); 벡터 d의도함수는각구간의중점에서의값이므로 >> xm=(x(1:n-1)+x(2:n))./2; 중점에대한 x 값을구한다. >> xa=0:.01 :.8 ; >> ya=25-400*xa+3*675*xa.^2-4*900*xa.^3+5*400*xa.^4; 도함수의해석값을구한다. >> xplot(xm, d, 'o', xa, ya) 수치해와해석해를그린다.
예제 21.4 (diff 을사용한미분 ) (4/4) 정확한도함수값 ( 실선 ) 과결과 ( 원 ) 의비교. MATLAB 의 diff 함수를사용해서계산한수치
예제 21.5 (gradient 를사용한미분 ) (1/2) Q. MATLAB 함수 gradient 를사용하여예제 21.4 의함수 를미분하라. Sol) >> f = @(x) 0.2+25*x-200*x.^2+675*x.^3-900*x.^4+400*x.^5; >> x = 0 : 0.1 : 0.8 ; >> y = f(x) ; >> dy=gradient(y, 0.1) gradient 함수를이용하여도함수를구한다. dy = Columns 1 through 5 10.8900 5.4400 1.5900 5.8400 8.5900 Columns 6 through 8 5.0400-4.8100-16.1600-21.3100 >> xa=0:.01 :.8 ; 수치해와해석해를그린다. >> ya=25-400*xa+3*675*xa.^2-4*900*xa.^3+5*400*xa.^4; >> xplot(x, dy, 'o', xa, ya)
예제 21.4 (diff 을사용한미분 ) (2/2) 결과는예제 21.4에서 diff 함수를사용해서얻은것보다정확하지못하다. 이는 diff(0.1) 을사용할때에비해두배가큰구간에서 gradient(0.2) 가실행되기때문이다. 정확한도함수값 ( 실선 ) 과 MATLAB의 gradient 함수를사용해서계산한수치결과 ( 원 ) 의비교.