Runge-Kutta 적분기
void rk4(double* y, double x,
void (*f)(double* dy, double* y, double x),
double h, int n);
// 운동방정식에 대한 연립 1차 미분방정식
void derv4(double* dz, double* z,double t);
#define DIM 4
void main()
{
int n;
double *y, h, t;
FILE *fp;
fp=fopen("problem8_2.dat", "w");
// 도함수 및 변수를 저장할 벡터를 생성한다.
// h : 적분기의 스탭크기를 0.01로 설정, n : 벡터 크기를 4로 설정
y = vector(DIM);
h = 0.01;
n = DIM;
t = 0.;
// 문제에서 주어진 값으로 초기화한다.
// y[0]=X, y[1]=V, y[2]=X, y[3]=V이다.
y[0] = 0.;
y[1] = 1.;
y[2] = 0.;
y[3] = 0.;
// 10초까지 0.01초 간격으로 적분
while(t<=10)
{
printf("t=%lf\n", t);
// rk4(Runge-Kutta함수)를 사용해 x과 x의 위치와 속도를 시간에 따라 파일에 저장
fprintf(fp, "%10.5lf %15.10lf %15.10lf %15.10lf %15.10lf\n", t, y[0], y[1], y[2], y[3]);
rk4(y, t, derv4, h, n);
t += h;
}
fclose(fp);
vector_free(y);
}
// derv4 : 문제1번에서 도출한 연립 1차 미분방정식이다.
void derv4(double *dz, double *z, double t)
{
dz[0] = z[1];
dz[1] = (z[2]-3*z[0]) / 4;
dz[2] = z[3];
dz[3] = z[0]-z[2];
}
결론
결과 및 고찰
수치해석을 수강하면서 선형 방정식 문제에 대해서 어떻게 하면 잘 표현을 할 수 있을까에 대해서 많은 내용을 배웠습니다. 이번 과제에서는 Runge-Kutta방법을 사용한 연립1차 미분방정식을 해결하는 문제를 풀어보았습니다. 이론시간에 배운 내용을 바탕으로 직접 C로 구현하는데 어려웠지만, 다행히 기본 함수를 위한 코드를 제공해주셔서 시간이 조금 덜 걸렸던 것 같습니다. Runge-Kutta를 이용해 풀어본 문제는 이론으로 배울 때는 어려운 점이 부분이 있었지만, 이를 C코드로 구현하는데 이론을 배울 때보다는 수월했던 것 같아 다행이었습니다. Runge-Kutta를 이용해 그래프를 도시하면서, 이전에 실습했던 수정오일러 방식과 비교해서 정상적으로 출력되어야하는 값과 좀 더 일치해지는 것을 확인할 수 있었습니다. 지금껏 배운 선형 방정식을 풀기위한 방법들을 사용해서 앞으로 전공으로 진행할 프로젝트에서도 같이 사용한다면 좋을 것 같다는 생각을 했습니다. 한 학기동안 수치해석이란 과목을 들으면서, 어려운 부분도 많이 있었지만, 그 만큼 얻어가는 지식도 많은 것 같습니다. 너무 수고하셨습니다.^_^