2007年10月31日 星期三

簡易用C語言求方程式的根

Standard
最近有不少人向我求救,原因不外乎是學校的報告、作業,而其中最令人覺得有趣的,就是關於數值方法的程式撰寫的問題。而講到用電腦處理數學運算,大多數人的直覺一定是使用 MATLAB 來實作其程式。是的,課堂上所教所講的就是 MATLAB,在處理各種數學問題和方程式,MATLAB 都有完整的函式指令和功能,只要熟悉,可以很輕易的能夠達成各種數學任務。不過最好笑的是,教授雖然教的是 MATLAB,但考的卻是 C 和 Fortran,要求學生的報告、作業也都是 C 和 Fortran,這讓一群學生都傻了眼,不知道該怎麼辦。

先暫且不管課本中是否附有 MATLAB 的程式碼,對於一群只上過幾門程式語言課的人,仍然無法輕易將程式轉用 C 或 Fortran 來寫。像是這樣一個用C語言求方程式的根的程式,就是許多人的問題:
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <conio.h>

#define F(x) sin(10*x)+cos(3*x)

double f(double x) {
return F(x);
}

int sign(double x) {
if (x<0)
return -1;
else if (x>0)
return 1;
else
return 0;
}

void incsearch(double xmin, double xmax, int ns)
{
int i,n;
double space,current,ans;

current = xmin;
space = (xmax - xmin) / (ns-1);

for (i=0;i<ns-1;i++) {
if (sign(f(current))!=sign(f(current+space))) {
printf("%.4f\t%.4f\n", current, current+space);
}
current += space;
}
}

int main(int argc, char *argv[])
{
int n;

printf("Number of subintervals: ");
scanf("%d", &n);
incsearch(3, 6, n);
printf("Orz...");

getch();
return 0;
}



這程式是求方程式:『f(x)=sin(10x)+cos(3x)』,x 在 3 到 6 之間的所有根。原理是將 x 軸做 n 次切割,並一一代入方程式以檢查根所在的範圍,所以 n 設越大,相對的答案範圍越精準。這部份我寫了個函式 incsearch() 尋找根,其利用正負來判斷根是否存在某區間內。

另外,本程式中實作了 MATLAB 中的 sign(),用來檢查數值是大於、等於或小於零。

後記

事實上求根要求的精準,還有更多方法可以使用,本文講的只是一種簡易求根的方法,想要知道更多,可以參考數值分析處理方面的書籍。亦或是日後又有人求救相關問題,在我幫忙完後有空也會在這記錄存檔起來。