Fred's blog

業界上常用的 Matlab 軟體，因為本身就是為了做數學運算而設計，所以在數值分析和處理方面非常快速好用，但是在 C 語言就中沒有那些好用的指令和數學函式庫了，得自己打造所有的運算流程和架構。舉例來說，使用 C 語言來做矩陣運算，程式開發者都是先『指定行列數量』建立一個陣列，再『指定行列數量』去對內部的值，一一做處理或高斯消去解方程式，這比起 Matlab 中用兩行就可以搞定一切運算，可說易用性相差甚遠。 所以就想到，要是我們能設計出一套好用的 C 語言矩陣函式庫，是否就可以省下這些重覆設計矩陣運算流程的時間？於是可以做一些實作，來完成這個任務。 基本的矩陣資料結構： typedef struct { int rows; int cols; int total; double *table; } Matrix; 其中定義了幾個基本的 functions 做矩陣的處理： 建立矩陣： Matrix *matrix_create(int rows, int cols) 矩陣最佳化： Matrix *matrix_optimize(Matrix *matrix) 解矩陣的方程式： Matrix *matrix_solution(Matrix *matrix) 顯示並列出矩陣內容： void matrix_print(Matrix *matrix) 清空矩陣： void matrix_clear(Matrix *matrix) 函式的程式碼： void matrix_clear(Matrix *matrix) { int i, j; double *ptr; /* initializing matrix */ ptr = matrix->table; for (i=0;i<matrix->total;i++) *(ptr++) = 0; } Matrix *matrix_create(int rows, int cols) { Matrix *matrix; /* allocate */ matrix = (Matrix *)malloc(sizeof(Matrix)); matrix->rows = rows; matrix->cols = col...

『梯型法則 Trapezoidal Rule』是積分方法的基本定理之一，從高中到大學的數學課程中都有提到，碰到工程這種需要大量求值的領域，我們就必需將這些積分方法寫成程式。 詳細程式碼： #include <stdio.h> #include <math.h> double my_eq_orig(double x) { return 0.2 + 25*x - 200*pow(x, 2) + 675*pow(x, 3) - 900*pow(x, 4) + 400*pow(x, 5); } double my_eq(double x) { return -400 + 4050*x - 10800*pow(x, 2) + 8000*pow(x,3); } double integrate_trapezoidal(double (*fx)(double), int n, double a, double b) { int i; double h; double x = a; double s = 0; h = (b - a) / n; for (i=1;i<=n-1;i++) { x += h; s = s + fx(x); } return h * ((fx(a) + fx(b)) / 2 + s); } void main(void) { /* real */ double real = integrate_trapezoidal(my_eq_orig, 100, 0, 0.8); /* n = 1 */ double real_i = integrate_trapezoidal(my_eq_orig, 1, 0, 0.8); /* average */ double fa = integrate_trapezoidal(my_eq, 100, 0, 0.8) / 0.8; printf("real value:\n%lf\n\n", real); printf("I:\n%lf\n\n", real_i); print...

依朋友們的要求，我實作了黃金比例搜尋法(Golden-Section Search) 的程式，用來尋找一個 f(x) 的最小值。這程式非常的簡短，是短短時間內寫出來的，所以並沒有考慮到效能、流程優化的部份，就像是 x 在某些狀況可以直接沿用到下次的 x1 或 x2，不需要將 x1、 x2 完全重新計算的這部份也是值得去探討的。 原始程式碼： #include <stdio.h> #include <math.h> #define F(x) ((pow(x,2)/10)-2*sin(x)) void main() { int i, n; double range_min = 0; double range_max = 4; double d, x1, x2; printf("Number of lines:"); scanf("%d", &n); printf("i xl f(xl) x2 f(x2) x1 f(x1) xu f(xu) d\n"); for (i=1;i<=n;i++) { d = 0.61803*(range_max - range_min); x1 = range_min + d; x2 = range_max - d; printf("%d %6.4f %6.4f %6.4f %6.4f %6.4f %6.4f %6.4f %6.4f %6.4f\n", i, range_min, F(range_min), x2, F(x2), x1, F(x1), ra...

最近有不少人向我求救，原因不外乎是學校的報告、作業，而其中最令人覺得有趣的，就是關於數值方法的程式撰寫的問題。而講到用電腦處理數學運算，大多數人的直覺一定是使用 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...