#include <conio.h> #include <cstdarg> #include <cstdio> #include "global.h" #define _STACK_CHECK_BOUNDARY_ template <class T, int N=512> class Stack { T& err (int n) { static char* s[] = { "Overflow","Underflow", "Out of range","Empty" }; puts (s[n]); getch(); return d[0]; } public: T d[N]; //資料項 int n; //元素個數 Stack() { reset(); } void reset() { n = 0; } #ifdef _STACK_CHECK_BOUNDARY_ void push (T x) { if (n<N) d[n++] = x; else err(0);} T pop() { return n? d[--n]: err(1); } T& operator*() { return n>0? d[n-1]: err(3); } T& operator[] (int i) { return 0<=i&&i<n? d[i]: err(2); } #else void push (T x) { d[n++] = x; } T pop() { return n? d[--n]: *d; } T& operator*() { return d[n-1]; } T& operator[] (int i) { return d[i]; } #endif }; //------------------------------------------------------------ // 使用方式: // // 1. 先 reset // 2. 以 add_vertex 依序由 0~N 號配置頂點 // 3. 用 add_edge 加入目前點的外接點群 // 4. 重複 2 3 步直到配置完畢 // //------------------------------------------------------------ struct Graph //圖形結構 { enum { MAX_V = 256, //節點最大個數 MAX_ES = 4096 //Edge Stack 的容量 }; struct Vertex { int* edge; //鄰接點群 int nEdge; //鄰接邊數 (接出去的邊數) }; Vertex V [MAX_V]; //Vertex stack int nV; //目前的頂點數 int curV; //目前正在做設置的頂點編號 int ES[MAX_ES]; //Edge Stack int nES; //目前的邊線數 void reset () { nV = nES = 0; } void add_vertex (int id) //id 為欲設置的頂點編號 { curV = id; if (id != nV) ERROR_MSG3_("頂點設置錯誤: %d 應為 %d", id, nV); V[curV].edge = &ES[nES]; V[curV].nEdge = 0; nV++; } void add_edge (int id) //id 為欲外接的頂點編號 { V[curV].nEdge++; ES[nES++] = id; } void add_edges (int v, int n...) { add_vertex (v); va_list p; va_start (p, n); while (n-- > 0) add_edge (va_arg (p, int)); } void print_edges () { int i, j; for (i=0; i<nV; ++i) { printf ("\n %2d: ", i); for (j=0; j<V[i].nEdge; ++j) printf ("%2d ", V[i].edge[j]); } } }; //------------------------------------------------------------ // 閉路搜尋: // // 使用 Tarjan 算法尋找 strongly connected nodes // 並略去元素個數為 1 或 2 的解 // //------------------------------------------------------------ class FindLoops //閉路搜尋算法 { typedef Stack<int,Graph::MAX_ES> TStack; Graph::Vertex* V; TStack stack; //保存正被遍歷的頂點群 int visit [Graph::MAX_V]; //存放頂點的 DFS 拜訪次序 bool bInStack [Graph::MAX_V]; //標示頂點是否在 stack 中 int order; //目前的拜訪次序 void tarjan (int v) //傳入起始節點編號 { static int cur; int min = visit[v] = order++; //min 存放頂點可追溯到的 //最前頂點的拜訪次序 int t; int n = V[v].nEdge; int* e = V[v].edge; stack.push (v); bInStack[v] = true; while (n-- > 0) { t = *e++; //取出鄰接點編號 if (visit[t] == -1) //以 DFS 拜訪鄰接點 tarjan (t); else cur = visit[t]; if (bInStack[t] && cur < min) min = cur; } if (visit[v] == min) { puts (""); do printf ("%2d ", t=stack.pop()), bInStack[t] = false; while (t != v); } } public: void operator() (Graph& graph) { stack.reset(); V = graph.V; order = 0; //-1 = 未被拜訪 memset (visit, -1, sizeof visit); for (int i=0; i<graph.nV; ++i) if (visit[i] < 0) tarjan (i); } };
2009年12月27日 星期日
閉路搜尋
2009年12月26日 星期六
前中後序運算式-2
以下程式說明了當運算子很多時，
最好用 reverse polish 的 stack 解法，而不要做 recursive descent，
否則 60 幾行可完成的 code 會變成好幾百行...
#include <iostream> #include <ctype.h> using namespace std;template <class T, int N=512> class Stack { T err (int n) { static char* s[] = { "Overflow","Underflow", "Out of range","Empty" }; puts (s[n]); getch(); return d[0]; } public: T d[N]; //資料項 int n; //元素個數 Stack() { n = 0; } void push (T x) { if (n<N) d[n++] = x; else err(0);} T pop() { return n? d[--n]: err(1); } T& operator*() { return n>0? d[n-1]: err(3); } T& operator[] (int i) { return 0<=i&&i<n? d[i]: err(2); } };class Expression { char *s, c, *p, buf[512]; void next () {c = *s++;} void put (char c1, char c2=0) { *p++ = c1; if (c2 != 0) *p++ = c2; } int expr () { if (!or()) return 0; while (c=='|') { next(); if(c != '|') return 0; next(); if(!or()) return 0; put ('|','|'); }return 1; } int or () { if (!and()) return 0; while (c=='&') { next(); if(c != '&') return 0; next(); if(!and()) return 0; put ('&','&'); }return 1; } int and () { if (!compar()) return 0; char op[3] = {0,0,0}; while (c=='<' || c=='>') { op[0] = c; next(); if(c == '=') {op[1]='='; next();} if(!compar()) return 0; put (op[0], op[1]); }return 1; } int compar () { if (!term()) return 0; for (char op; c=='+' || c=='-';) { op = c; next(); if(!term()) return 0; put (op); }return 1; } int term () { if (!fact()) return 0; for (char op; c=='*' || c=='/';) { op = c; next(); if (!fact()) return 0; put (op); }return 1; } int fact () { if (!expo()) return 0; while (c=='^') { next(); if (!expo()) return 0; put ('^'); }return 1; } int expo () { if (c == '!') { next(); if(!invers()) return 0; put ('!'); } else if (!invers()) return 0; return 1; } int invers () { if (isalpha(c)) { put (c); next(); return 1; } else if (c=='(') { next(); if (expr() && c == ')') { next(); return 1; } }return 0; } //--------------------------------- struct Node { char c[3]; Node *L,*R; Node () { c[1]=c[2]=0; L=R=0; } Node (char c1, char c2=0, Node* l=0, Node* r=0) { c[0]=c1; c[1]=c2; c[2]=0; L=l; R=r; } Node (char c1, Node* l, Node* r) { c[0]=c1; c[1]=0; L=l; R=r; } }; Node* build_exptree (char* postfix) { Stack <Node*> s; Node *o, *L, *R; char c1, c2, *p = postfix; while(*p) { if(isalpha(*p) || isdigit(*p)) { s.push (new Node(*p)); } else { c1 = *p; c2 = 0; L = 0; R = s.pop(); if (c1 != '!') L = s.pop(); if (c1=='&'||c1=='|') c2 = *++p; else if ((c1=='<'||c1=='>')&& p[1]=='=') c2 = *++p; s.push (new Node (c1, c2, L, R)); } p++; } return s.pop(); } void preorder (Node *x) { if(x != NULL) { cout << x->c; preorder (x->L); preorder (x->R); } } void inorder (Node *x) { if(x != NULL) { preorder (x->L); cout << x->c; preorder (x->R); } } void postorder (char* exp) { p = buf; s = exp; next(); expr(); put(0); } //----------------------------------- public: void postfix (char* exp) { postorder (exp); cout << buf; } void prefix (char* exp) { postorder (exp); preorder (build_exptree (buf)); } }; void main (/* daviddr */) { char* s[] = { "(a-b*(c+d)^e)&&((f+g^h)-i*j+k/l)", "!((a^b-c/d+e/f))||(g-(h^i+j*k)/l+m)", "(a-b/c)*d^f<=(g+h*i)/j+k^l-m" }; Expression exp; for (int i=0; i<3; i++) { puts(""); exp.postfix(s[i]); puts(""); exp.prefix (s[i]); } }
2009年12月25日 星期五
有向鄰接圖資料結構
若頂點總數為 V
所有頂點的鄰接點數總和為 N
常見的結構有:
1. Adjacency Matrix
使用空間 V * V
2. Adjacency List:
記錄 out edge
若串列指標總和大小為 N
使用空間 V * (2*(1+N))
3. Adjacency Multilists:
用 2 個 list 記錄 out edge 與 in edge
使用空間 2 * V * (2*(1+N))
4. Sequential Representation 或 Index Table:
當頂點數多但 edge 少時，相當省空間與時間。
使用空間 V + (V+N)
前頭的 V 為索引表大小
所有頂點的鄰接點數總和為 N
常見的結構有:
1. Adjacency Matrix
使用空間 V * V
2. Adjacency List:
記錄 out edge
若串列指標總和大小為 N
使用空間 V * (2*(1+N))
3. Adjacency Multilists:
用 2 個 list 記錄 out edge 與 in edge
使用空間 2 * V * (2*(1+N))
4. Sequential Representation 或 Index Table:
當頂點數多但 edge 少時，相當省空間與時間。
使用空間 V + (V+N)
前頭的 V 為索引表大小
2009年12月24日 星期四
高維下的 kd-tree
http://www.codeproject.com/KB/architecture/KDTree.aspx
Anton Milev 文中寫到：
For small dimensions the KD-tree can be much faster then sequential search, however for 10-dimensional space the checked nodes for rises to about 25% for. It is possible to reduce this number 2-3 times if the tree is balanced but the tendency of exponential rise of the checked nodes makes the KD unpracticle above 15 dimensions. The second curve shows how KD / Brute ratio changes with space dimension, Brute stands for sequential search, it checks all nodes. Practiclly, for 10000 points, the KD tree becomes as fast as the sequential search for dimensions bigger then 10.
With the analyses made in this article I showed that the standard KD-tree is not good for dimensions beyond 10. In my next article about KD trees and the nearest neighbours search in N-dimensions I will introduce an extention of the quad tree for N-dimensional data - QD-tree, it is faster then KD and allows run-time balancing. However for big dimensions it again becomes slow. This problem is well known in the computational geometry, this would be so for any data structure based only on tree. From another hand it is not possible to make a multidimensional linked-list and to include it in the tree. The problem is quite complicated and requires new kinds data structures.
2009年12月23日 星期三
Voronoi Diagram
Voronoi Diagram / Delaunay Triangulation
http://www.cs.cornell.edu/Info/People/chew/Delaunay.html
+Crust
http://www.cs.unc.edu/~snoeyink/demos/crust/home.html
http://www.cs.cornell.edu/Info/People/chew/Delaunay.html
+Crust
http://www.cs.unc.edu/~snoeyink/demos/crust/home.html
2009年12月18日 星期五
字串分行且左右對齊
為使編排較美觀，我根據單字長度，由短而長依次調整字距：
//以 a 為 primal key, b 跟著 a 一塊兒做 sorting void sort (int* a, int* b, int n) { int t, u, i, j, g; for (g=n; g>0; g= (g>1 && g<5)? 1: 5*g/11) for (i=g; i<n; ++i){ t = a[i]; u = b[i]; for (j=i; j>=g && t<a[j-g]; j-=g) { a[j] = a[j-g]; b[j] = b[j-g]; } a[j] = t; b[j] = u; } } int main (/* -.. .- ...- .. -.. -.. .-. */) { int index[512]; //反向索引 int wlen [512]; //單詞長度 int stuff[512]; //填充空白數 int i = 0, nStuff; //i 為累計單詞數 int j, n, len = 0; int width = 40; //width = 指定欄寬 //欄寬需大於最大單詞長度! char in[] = "In computing, C is a general-purpose, block structured," " procedural, imperative computer programming language " "developed in 1972 by Dennis Ritchie at the Bell Telephone" "Laboratories for use with the Unix operating system."; char *s = in, *p = in; char out[4096], *o = out; while (*p) { stuff[i] = 0; index[i] = i; for (n=0; *p && ' '!=*p; ++n, ++p); //n = 單詞長度 if (' '==*p) { ++p; ++n; } len += (wlen[i] = n); //累計列長 if (len > width) { nStuff = width - len + n; //算總填充數 sort (wlen, index, i-1); //根據字長排序 j = 0; while (nStuff--) { //算個別填充數 ++stuff[j]; if (++j >= i-1) j = 0; } sort (index, stuff, i-1); //將 stuff 照 //原始字序排列 for (j=0; j<i-1; ++j) { while (' '!=*s) *o++ = *s++; do *o++ = ' '; while (stuff[j]--); ++s; } while (*s && ' '!=*s) *o++ = *s++; if (' ' == *s) { *o++ = '\n'; ++s; } index[0] = stuff[0] = 0; wlen[0] = len = n; i = 0; } ++i; } for (p-=len; *p; *o++ = *p++); *o++ = *p++; puts (out); return system ("pause"); }
寬度為 40 時，輸出:
In computing, C is a general-purpose, block structured, procedural, imperative computer programming language developed in 1972 by Dennis Ritchie at the Bell Telephone Laboratories for use with the Unix operating system.
void main
很久以前寫的..
久到幾乎忘了曾寫過一堆 void main..
void main() 為編譯器廠商自訂的簡便用法，
在不同 OS 與 compiler 上不具通用性。
許多早期的 C Compiler，無論是 void main()、int main()、或 main()
最後在組語符號表中的都是 _main，於是會發現在某些編譯器 (如 K&R) 上，
void main(), double main() 或是 int main[]= {0xcc} 也能順利編譯 :-)
( ) 中不加參數宣告時.. 多數編譯器視為 (void)，
亦有視為不定參數者，不過那是陳年往事了，現在多不這麼做。
x86 系統上，最後執行 xor eax, eax，(接著 pop edi、esi、ebx，
回復 stack 後 ret) 將 0 放在 EAX 暫存器裡傳回給作業系統。
C89/C99 中規定 main 需指明傳回值，
因此新式、標準的編譯器不再支援 main () { ... } 的寫法。
main 可移植的正規寫法有 2 種:
int main (void);
int main (int argc, char *argv[]);
Platform-dependent 格式有 2 種:
(1) UNIX、DOS、Windows 可用:
int main(int argc, char *argv[], char *envp[]);
(2) Mac OS X 與 Darwin 可用:
int main(int argc, char *argv[], char *envp[], char *apple[]);
envp 中放的就是用 getenv 取出的環境變數。
各參數意義 google 上可查到許多介紹。
取得正在執行的程式列表
一年多前的簡答，或許會需要用到，記錄一下。
#include <windows.h> #include <tlhelp32.h> using namespace System; int main (array<String^>^) { HANDLE hProcs = CreateToolhelp32Snapshot (TH32CS_SNAPPROCESS, 0); if (!hProcs) return 0; PROCESSENTRY32 pe = {0}; pe.dwSize = sizeof (PROCESSENTRY32); if (Process32First (hProcs, &pe)) do Console::WriteLine (L"PID: {0}\nName: {1}\n", pe.th32ProcessID, gcnew String (pe.szExeFile)); while (Process32Next (hProcs, &pe)); Console::ReadKey(); return 0; }
2009年12月17日 星期四
JNI
使用 JNI 最保險的方式是透過 LoadLibrary 導入 jvm.dll。
首先設好 include 和 lib path，這個網路上資料很多，可自行查找。
錯誤排解:
[1] 發生 Java 異常 EXCEPTION_ACCESS_VIOLATION (0xc0000005):
表示 jvm.dll 位置有誤。
直接指定到 /jre/bin/client/ 去 LoadLibrary 就行，
任何「多餘」的拷貝、指定動作都可能引起此異常。
[2] JNI_CreateJavaVM 傳回 JNI_ERR:
「可能」是先前安裝過多個 JDK 或 JRE 版本，
因未移除，或未移除乾淨，互相干擾所致。
解法是直接 load jvm.dll。
[3] Dev C++ 如何加 -ljvm?
沒時間試驗，但至少 lib path 要設對。
Test.java
供 GetStaticMethodID 識別調用。
首先設好 include 和 lib path，這個網路上資料很多，可自行查找。
錯誤排解:
[1] 發生 Java 異常 EXCEPTION_ACCESS_VIOLATION (0xc0000005):
表示 jvm.dll 位置有誤。
直接指定到 /jre/bin/client/ 去 LoadLibrary 就行，
任何「多餘」的拷貝、指定動作都可能引起此異常。
[2] JNI_CreateJavaVM 傳回 JNI_ERR:
「可能」是先前安裝過多個 JDK 或 JRE 版本，
因未移除，或未移除乾淨，互相干擾所致。
解法是直接 load jvm.dll。
[3] Dev C++ 如何加 -ljvm?
沒時間試驗，但至少 lib path 要設對。
Test.java
public class Test { public static void hello () { System.out.print ("hello~"); } public static void main (String[]s) {} }以 javap -s -p Test 觀察各 method 的 Signature，
供 GetStaticMethodID 識別調用。
#undef UNICODE #pragma comment (lib,"jvm.lib") #include "windows.h" #include "jni.h" #include <iostream> using namespace std; #define CHECK(o) if(o); else{\ cout<<__FILE__<<" ("<<__LINE__<<"): not satisfy "<<#o ;\ exit (getchar());\ } typedef jint (WINAPI* PFunCreateJavaVM) (JavaVM**, void**, void*); int main () { JavaVMOption opt[4] = { {"-Djava.compiler=NONE",0}, {"-Djava.class.path=." ,0}, }; JavaVMInitArgs vm_args = { JNI_VERSION_1_2, 2, opt, JNI_TRUE //opt = JRE 運行參數 }; JNIEnv* env = 0; JavaVM* jvm = 0; HINSTANCE hInst = LoadLibrary ("C:/Java/jre/bin/client/jvm.dll"); CHECK (hInst != 0); PFunCreateJavaVM CreateJavaVM = (PFunCreateJavaVM) GetProcAddress (hInst, "JNI_CreateJavaVM"); CreateJavaVM (&jvm, (void**)&env, &vm_args); CHECK (env != 0); jclass cls = env->FindClass ("Test"); //查找 Java 類別 CHECK (cls != 0); jmethodID pfn = env->GetStaticMethodID (cls, "hello", "()V"); CHECK (pfn != 0); env->CallStaticVoidMethod (cls, pfn); if (env->ExceptionOccurred ()) { env->ExceptionDescribe (); } jvm->DetachCurrentThread (); jvm->DestroyJavaVM (); return getchar(); }
2009年12月12日 星期六
堆疊式線段追蹤
用 recursion 在線段複雜時會造成 stack overflow，
因此需要 heap 解。
LineTrace 提供了 singleton 介面，但不建議在大型專案中使用，
因為它除了有 global variable 的缺點外，
和其他類別合作時，會引入建構解構上的相依性，難以 debug。
應盡量把 LineTrace 設計成 class component 來用。
#include "global.h" #include "help_fn.h" #include "new_arr.h" //------------------------------------------------------------- // state 紀錄圖中每一點的追蹤狀態: // // ６７８ // ５米１ // ４３２ // // 狀態 = 0 表示非邊線點 //------------------------------------------------------------- class LineTrace { enum {W=512, H=512}; //狀態圖的最大寬高 int w, h; //檢測邊界 IplImage* out; //輸入邊線圖, 同時用來存放輸出結果 IplImage* mark; //阻隔點, 遇阻隔點則不再追蹤下去 int nPos; //堆疊頂端指標 int* posx; //紀錄追蹤位置的堆疊 int* posy; // uchar* state; //標示圖中每一點的追蹤狀態 static LineTrace* self; //Singleton void check_vaild() //確認資料的有效性 { if (out == 0) ERROR_MSG_("out == 0"); if (mark== 0) ERROR_MSG_("mark == 0"); if (out->width > W) ERROR_MSG3_("影像寬度: %d > %d", out->width, W); if (out->height > H) ERROR_MSG3_("影像高度: %d > %d", out->height, H); } public: static LineTrace& instance() { if (self == 0) self = new LineTrace; return *self; } static void release() { delete self; } LineTrace() { posx = new int [W*H]; posy = new int [W*H]; state = new uchar[W*H]; } ~LineTrace() { delete[] posx; delete[] posy; delete[] state; } void setImage (IplImage* out_, IplImage* mark_) { out = out_; mark = mark_; check_vaild(); } void trace (int ox, int oy) // (ox,oy) = trace 的起始座標 { static CvScalar red = cvScalar (0,255,0); int x, y; //目前 trace 到的點座標 int dir; //目前要 trace 的方位 CvScalar c, c2; check_vaild(); nPos = 0; #define PUSH(px,py) \ nPos++;\ posx[nPos] = x = px;\ posy[nPos] = y = py;\ dir = state[nPos] = 1; #define POP \ nPos--;\ x = posx[nPos];\ y = posy[nPos];\ dir = state[nPos]; PUSH (ox,oy); //邊界衛兵 :-) PUSH (ox,oy); //將起始座標推入堆疊 //開始進行線段追蹤 while (nPos > 1) { if (nPos >= W*H) ERROR_MSG_("stack overflow!"); if (dir == 1) { //若是第一次檢測某個座標點 c = cvGet2D (out, y, x); //就取出顏色值 c2 = cvGet2D (mark, y, x); if (c.val[0] == 0 || //若不是邊緣線 c2.val[0] > 0){ //或是遇到阻礙點 POP; //便回復先前的狀態 continue; //退回上個節點 } cvSet2D (out, y, x, red); //標示追蹤過的線 } if (dir <= 0) { //這是不可能發生的 puts ("error!"); getch(); } if (dir > 8) { //若所有方向都檢測過了 POP //便回復先前的狀態 continue; //退回上個節點 } else { state[nPos]++; //遞增目前節點的方向值 switch (dir) { //移至特定方位的鄰居點 case 1: x++; break; case 2: x++; y++; break; case 3: y++; break; case 4: x--; y++; break; case 5: x--; ; break; case 6: x--; y--; break; case 7: y--; break; case 8: x++; y--; break; } PUSH (x,y); //將鄰居點推入堆疊頂端 } //進行邊線追蹤 } #undef PUSH #undef POP } }; LineTrace* LineTrace::self = 0; #define ef else if //------------------------------------------------------------- // 說明: // 壓著左鍵移動滑鼠 - 畫出線條 // 按中間的滾輪鍵 - 加入阻隔點 // 在線段上按右鍵 - 進行線段追蹤 // 按空白鍵 - 重置畫面 // 按 ESC - 離開程式 //------------------------------------------------------------- struct Painter { static LineTrace tracer; //線段追蹤器 static IplImage* src; //畫板 static IplImage* mark; //阻隔點 static CvPoint pt; //先前的點座標 Painter (char* path = 0) //可傳入黑白圖檔 { if (path) src = cvLoadImage (path, CV_LOAD_IMAGE_GRAYSCALE); else src = cvCreateImage (cvSize (400,400), 8, 3); mark = cvCreateImage (cvGetSize (src), 8, 3); reset(); tracer.setImage (src, mark); cvNamedWindow ("image", 1); cvSetMouseCallback ("image", on_mouse, 0); } ~Painter() { //釋放資源 cvDestroyWindow ("image"); cvReleaseImage (&src); cvReleaseImage (&mark); } static void reset() //清除畫面 { cvZero (src); cvZero (mark); } static void on_mouse (int event, int x, int y, int flags, void*) { static CvScalar white = CV_RGB (255,255,255); if (src == 0) return; if (event == CV_EVENT_LBUTTONUP) //如果放開了滑鼠左鍵 //!(flags & CV_EVENT_FLAG_LBUTTON)) //或按下了其他鍵 pt = cvPoint (-1,-1); //便重置描繪起始點 ef (event == CV_EVENT_MBUTTONUP) { cvCircle (mark, cvPoint(x,y), 4, white, -1); cvCircle (src, cvPoint(x,y), 4, CV_RGB(0,180,255), -1); show(); } ef (event == CV_EVENT_RBUTTONUP) { puts ("start tracing"); for (int j,i=-1; i<2; ++i) for (j=-1; j<2; ++j) tracer.trace (x+i, y+j); puts ("end tracing"); show(); } ef (event == CV_EVENT_LBUTTONDOWN) pt = cvPoint (x,y); ef (event == CV_EVENT_MOUSEMOVE && //若移動滑鼠 (flags & CV_EVENT_FLAG_LBUTTON) ) //且壓著左鍵 { CvPoint p = cvPoint (x,y); if (pt.x < 0) pt = p; cvLine (src, pt, p, white, 1, 8, 0); pt = p; show(); } } static void show() { cvShowImage ("image", src); } }; LineTrace Painter::tracer; IplImage* Painter::src; IplImage* Painter::mark; CvPoint Painter::pt; int main() { Painter painter; painter.show(); while (true) { int key = cvWaitKey(0); if (key == 27) break; ef (key == ' ') { painter.reset(); painter.show(); } } }
執行結果：
畫了張圖做測試~~~
2009年12月10日 星期四
老鼠走迷宮
步驟：
[1] 往上下左右拜訪鄰接點。
[2] 若鄰點可通行，便遞迴至 [1]。
內有暴力 rand() ...
#include <windows.h> #include <ctime> #include <iostream> using namespace std; const int W = 23, H = 32; //迷宮長寬 wchar_t m[W+1][H+1] = {}; //迷宮本體 int X, dX[] = {0,-2,0,2}, d; //方位角 int Y, dY[] = {2,0,-2,0}; bool bEnd = false; //是否已走到終點 void make (int x, int y) { m[x][y] = L' '; for (int n = 16; n--; ) { d = rand() % 4; X = x + dX[d]; Y = y + dY[d]; if (0<X && X<W && 0<Y && Y<H && m[X][Y] == L'■') m[x+dX[d]/2][y+dY[d]/2] = L' ', make (X,Y); } } void show_step (int x, int y, wchar_t ch) { COORD c = {y*2, x}; //顯示過程 SetConsoleCursorPosition ( GetStdHandle (STD_OUTPUT_HANDLE), c); wcout << (m[x][y] = ch); Sleep (50); } void visit (int x, int y) //尋路 { if (x==X && y==Y) bEnd = true; //已抵達終點 else if (m[x][y] == L' ' && !bEnd) { show_step (x, y, L'˙'); for (int i=0; i<4; ++i) visit (x+dX[i]/2, y+dY[i]/2); if (!bEnd) show_step (x, y, L' '); } } int main (/* daviddr 091211 */) { srand (time(0)); wcout.imbue (locale ("cht")); for (d=W; d--;) wmemset (m[d], L'■', H)[H] = L'\n'; make (1,1); X = W-2; //記錄終點座標 Y = H-1; m[1][0] = m[X][Y] = L'◎'; wcout << *m; visit (1,1); getchar(); }執行結果：
■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■ ◎˙˙˙˙˙■˙˙˙˙˙■ ˙˙˙˙˙˙˙■ ■ ■■■■■˙■˙■■■˙■ ■■■■■˙■■■■■˙■ ■■■■■■■ ■ ■˙˙˙˙˙■˙˙˙■˙■ ■˙■ ■˙■ ■ ■˙˙˙˙˙■ ■˙■■■■■ ■˙■˙■■■■■■■˙■■■ ■˙■ ■ ■˙■■■˙■ ■˙˙˙■ ■˙■˙■˙˙˙˙˙■˙■ ■˙■ ■˙■ ■˙■ ■■■˙■ ■■■˙■˙■˙■■■˙■˙■ ■■■˙■■■ ■˙■ ■˙■ ■˙˙˙■ ■˙■˙˙˙■ ■˙■˙■ ■˙˙˙■ ■˙■ ■˙■ ■˙■■■■■■■˙■■■■■ ■˙■˙■ ■˙■■■■■■■˙■ ■˙■ ■˙˙˙˙˙˙˙˙˙■ ■˙˙˙■˙˙˙■ ˙˙˙˙˙■˙˙˙■ ■■■■■■■■■■■■■ ■ ■■■■■˙■■■■■˙■■■■■˙■■■ ■ ■ ■ ■˙˙˙˙˙˙˙■˙˙˙˙˙■ ■ ■ ■■■■■ ■ ■■■■■■■■■ ■■■■■■■■■˙■■■■■ ■ ■ ■ ■˙˙˙■ ■˙˙˙ ■˙■ ■ ■ ■■■ ■■■ ■˙■˙■■■■■■■ ■˙■˙■ ■˙■■■■■ ■ ■ ■ ■ ■ ■˙■˙˙˙˙˙˙˙■ ■˙■˙■ ■˙˙˙˙˙■ ■ ■ ■ ■■■ ■■■˙■■■■■■■˙■■■˙■˙■■■■■■■˙■ ■ ■ ■ ■ ■˙˙˙˙˙■ ■˙˙˙˙˙■˙˙˙˙˙˙˙˙˙■ ■ ■■■ ■ ■■■ ■■■■■˙■ ■■■■■■■■■■■■■■■■■ ■ ■ ■ ■ ■ ■˙˙˙■ ˙˙˙■ ˙˙˙■˙˙˙˙˙■ ■ ■ ■ ■ ■■■ ■˙■■■■■■■˙■˙■■■˙■˙■˙■■■˙■ ■ ■ ■ ˙˙˙˙˙˙˙˙˙■˙˙˙˙˙■˙˙˙■ ˙◎ ■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■■
使用 shuffle 可以隨機拜訪鄰接點，
而 make 中的 for 只要施行 4 次 (轉折點也會變多...(吧?))。
隨機拜訪鄰點，或許能稍微縮短抵達終點的平均時間
(但也可能更久，需視迷宮結構而定)。
#include <algorithm> void visit (int x, int y) //尋路 { if (x==X && y==Y) { //已抵達終點 bEnd = true; return; } if (m[x][y] == L' ' && !bEnd) { show_step (x, y, L'˙'); int i, D[] = {0,1,2,3}; //鄰點序列 random_shuffle (D, D+3); //隨機拜訪鄰點 for (i=0; i<4; ++i) visit (x+dX[D[i]]/2, y+dY[D[i]]/2); if (!bEnd) show_step (x, y, L' '); } }
OEM懷舊顯示 :-)
#include <windows.h> void setCodePage (int n...) { va_list p; va_start (p, n); while (!SetConsoleOutputCP(n)) n = va_arg (p, int); } void main () { int W=20, i; char *p, const*const a = "XXXXXXXXXXXXXXXXXXXX" "X X X" "X XXX X X XXX XX X X" "X X X X X X X" "X X X X X X X" "X X X X XXX X X" "X X X X X X X" "X X X X XXX X X X" "X X X" "X XXXXXXXXXXXXXXXXXX"; setCodePage (437, 855, 10000); for (p=a, i=0; *p; i++) { if (i>=W) { puts(""); i=0; } printf ("%c", *p++=='X'?219:' '); } }
2009年12月9日 星期三
無框視窗
使用 Dev C++ 的話，連結器欄需加入：-lgdi32
#include <windows.h> #include <iostream> using namespace std; extern "C" WINBASEAPI HWND WINAPI GetConsoleWindow(); void set_size (int w, int h) { int cy = GetSystemMetrics (SM_CYCAPTION), bx = GetSystemMetrics (SM_CXBORDER), by = GetSystemMetrics (SM_CYBORDER), fx = GetSystemMetrics (SM_CXFIXEDFRAME), fy = GetSystemMetrics (SM_CYFIXEDFRAME); if (w > 640) w = 640; HANDLE hOut = GetStdHandle (STD_OUTPUT_HANDLE); HWND hwnd = GetConsoleWindow (); SMALL_RECT rc = {0, 0, w/8 ,h/16}; COORD size = {w/8+1 ,h/16+1}; int dx = fx+bx; int dy = cy+fy+by; SetConsoleWindowInfo (hOut,1,&rc); //變更視窗大小 SetConsoleScreenBufferSize (hOut, size); SetWindowRgn (hwnd, CreateRectRgn (dx, dy, w-dx, h-dy), 1); } int main (/* 091210 */) { set_size (490, 690); cout << "frameless window"; return getchar(); }
Simple K-Means
參考了 Adrien Le Batteux 和 OpenCV 的 kmeans2 範例。
不過寫的落落長，實在不怎麼 simple @@
點集分群結果：
程式碼：
struct FontObject { char str[256]; CvFont font; //字型物件 CvPoint pos; CvScalar color; FontObject (double sw=0.25, double sh=0.25, int thickness = 1) { color = cvScalar (255,255,255); cvInitFont (&font, CV_FONT_HERSHEY_SIMPLEX, sw,sh, 0, thickness, CV_AA); } void show (IplImage* out, int x, int y) { pos.x = x; pos.y = y; sprintf (str, "(%d,%d)", x, y); cvPutText (out, str, pos, &font, color); } void show (IplImage* out, int x, int y, char* txt) { pos.x = x; pos.y = y; cvPutText (out, txt, pos, &font, color); } }; FontObject fontObj (0.75, 0.75, 2); template <typename T> void shuffle (T* a, int n) //簡易洗牌 { while (n > 1) { int j = rand() % n; n--; T t = a[j]; a[j] = a[n]; a[n] = t; } } template <typename T> inline void set_vector (T* out, T* in, int n) { while (n-- > 0) out[n] = in[n]; } template <typename T> inline void add_vector (T* out, T* in, int n) { while (n-- > 0) out[n] += in[n]; } template <typename T> inline void div_vector (T* v, int divisor, int n) { while (n-- > 0) v[n] /= divisor; } template <typename T> inline double cal_distance (T* a, T* b, int n) { double sum = 0, diff; while (n-- > 0) { diff = a[n] - b[n]; sum += diff * diff; } return sqrt (sum); } template <typename T = float> struct Kmeans { int nIteration; //最大迭代次數 int nData; //資料項數 int* label; //資料項的分類標籤 int* ccount; //群心涵蓋的資料個數 T** data; //指向資料集合 T** center; //群心集合 T* new_center; //計算暫用的群心 T* mean; //資料項的幾何中心 int nCenter; //群心個數 int dim; //資料維度 double error; //群心更新時的最大移動量 Kmeans() { srand (time(0)); label = ccount = 0; center = data = 0; new_center = mean = 0; } ~Kmeans() { release(); } private: void set (T* out, T* in) {set_vector (out, in, dim);} void add (T* out, T* in) {add_vector (out, in, dim);} void div (T* v, int divisor) {div_vector (v, divisor, dim);} double distance (T* a, T* b) {return cal_distance (a, b, dim);} void release() //釋放資源 { if (label) del_arr (label); label = 0; if (ccount) del_arr (ccount); ccount = 0; if (center) del_arr (center); center = 0; if (mean) del_arr (mean); mean = 0; if (new_center) del_arr (new_center); new_center = 0; } void init_centers () { if (nCenter < 1) ERROR_MSG_("群心數目需 > 1"); if (nCenter > nData) nCenter = nData; release(); //釋放資源 new_center = (T*) new_arr (sizeof(T), 1, dim); mean = (T*) new_arr (sizeof(T), 1, dim); center = (T**) new_arr (sizeof(T), 2, nCenter, dim); ccount = (int*) new_arr (sizeof(int), 1, nData); label = (int*) new_arr (sizeof(int), 1, nData); int i; for (i=0; i<nCenter; ++i) //隨機指定樣本點為群心 set (center[i], data [rand() % nData]); for (i=0; i<nData; ++i) //計算資料項的幾何中心 add (mean, data[i]); div (mean, nData); } void cluster_assignment () //將樣本分群 { int label_id, j, k; double min_dist, dist; memset (ccount, 0, nData* sizeof(int)); for (j=0; j<nData; ++j) { label_id = 0; min_dist = distance (center[0], data[j]); for (k=1; k<nCenter; ++k) { dist = distance (center[k], data[j]); if (dist < min_dist) { label_id = k; min_dist = dist; } } label[j] = label_id; //歸類 ccount[label_id]++; //統計群心涵蓋的資料個數 } } void update_cluster () //重新估算 centroids 位置 { int j, k, total; double dist = 0; error = 0; for (k=0; k<nCenter; ++k) { memset (new_center, 0, dim* sizeof(T)); total = 0; for (j=0; j<nData; ++j) if (label[j] == k) { add (new_center, data[j]); total++; } if (total > 0) { div (new_center, total); dist = distance (center[k], new_center); if (error < dist) error = dist; set (center[k], new_center); } else set (center[k], mean); } } public: //算法主流程: nCenter = 群心數目, // data = 資料集合 = data [nData_] [dim_] void do_cluster (int nCenter_, double error_threshold, T** data_, int nData_, int dim_) { nCenter = nCenter_; //設定群心數目 data = data_; //指向資料集合 nData = nData_; //設定資料個數 dim = dim_; //設定資料維度 init_centers (); //初始化群心位置 for (int i=0; i<nIteration && error > error_threshold; ++i) { cluster_assignment (); //將樣本分群 update_cluster (); //重新估算 centroids 位置 } } //單步執行 clustering, 回傳群心最大移動量 void setup (int nCenter_, T** data_, int nData_, int dim_) { nCenter = nCenter_; //設定群心數目 data = data_; //指向資料集合 nData = nData_; //設定資料個數 dim = dim_; //設定資料維度 init_centers (); //初始化群心位置 } double iterater_cluster () { cluster_assignment (); //將樣本分群 update_cluster (); //重新估算 centroids 位置 return error; } }; //----------------------------------------------------------------- template <class T> struct KMeanTest { enum {CD = 3, MAX_CLUSTERS = (CD+1)*(CD+1)*(CD+1), //最大群心數 W=700, H=700}; //展示影像的寬高 bool bEnd; CvScalar color_tbl[MAX_CLUSTERS]; IplImage* img; //2D展示影像 Kmeans<T> kmean; T** data; //資料集合 int nCenter; //群心個數 int nData; //資料項數 int dim; //資料維度 private: void release() //釋放資源 { if (data) del_arr (data); data = 0; } void generate_color_table() { int r, g, b, i=0, d = 255/CD; for (r=0; r<256; r+=d) for (g=0; g<256; g+=d) for (b=0; b<256; b+=d) color_tbl[i++] = CV_RGB (r,g,b); shuffle (color_tbl+1, MAX_CLUSTERS-2); } void generate_sample () { release(); data = (T**) new_arr (sizeof(T), 2, nData, dim); for (int i=0; i<nData; ++i) { data[i][0] = rand() % W; data[i][1] = rand() % H; } } public: KMeanTest() { generate_color_table(); cvNamedWindow ("clusters", 1); data = 0; img = cvCreateImage (cvSize (W, H), 8, 3); } ~KMeanTest() { cvDestroyWindow ("clusters"); cvReleaseImage (&img); release(); } void show() { static CvScalar white = CV_RGB (255,255,255); static char txt[9]; int i, key; cvZero (img); for (i=0; i<nData; ++i) //畫出資料項 { CvPoint pt = cvPoint (data[i][0], data[i][1]); CvScalar color = color_tbl [kmean.label[i] + 1]; cvCircle (img, pt, 2, color, CV_FILLED, CV_AA, 0); } for (i=0; i<nCenter; ++i) //畫出群心位置 { T* c = kmean.center[i]; cvCircle (img, cvPoint (c[0],c[1]), 3, white, CV_FILLED, CV_AA, 0); sprintf (txt, "%d", i); //顯示群心標號 fontObj.show (img, c[0]+3, c[1], txt); } cvShowImage ("clusters", img ); key = cvWaitKey(0); if (key == 27 || key == 'q' || key == 'Q') //ESC bEnd = true; } void start (int nCenter_, int nData_, int dim_, double threshold) { bEnd = false; nCenter = nCenter_; nData = nData_; dim = dim_; if (nCenter > MAX_CLUSTERS-2) { cout << "群心數目不得大於 " << MAX_CLUSTERS-2; return; } generate_sample (); kmean.setup (nCenter, data, nData, dim); while (!bEnd && kmean.iterater_cluster() > threshold) show(); } }; int main () { KMeanTest<double> test; test.start(24, 10000, 2, 5); }
2009年12月6日 星期日
simple iostream
改自：http://okmij.org/ftp/packages/iostream.h
#include <stdio.h> static class Iostream {} cin, cout; static class Endl{} endl; Iostream& operator << (Iostream& o, Endl& data) { puts (""); return o; } Iostream& operator << (Iostream& o, void* data) { printf ("%x", data); return o; } Iostream& operator << (Iostream& o, bool data) { printf ("%s", data? "true": "false"); return o; } Iostream& operator >> (Iostream& i, char* data) { scanf ("%s", data); return i; } #define DEF_OS(T,fmt) \ Iostream& operator << (Iostream& o, T data) {\ printf ("%"#fmt,data);\ return o;\ }\ Iostream& operator >> (Iostream& i, T& data) {\ scanf ("%"#fmt,&data);\ return i;\ } DEF_OS (const char* , s) DEF_OS (char* , s) DEF_OS (char , c) DEF_OS (int , d) DEF_OS (float , g) DEF_OS (double , lg) DEF_OS (long , ld) DEF_OS (long long , ll) DEF_OS (unsigned , u) DEF_OS (unsigned long, lu) int main () { char s[64]; double d; cin >> s >> d; cout << "input: " <<s <<' ' << d <<endl <<&d; }
2009年11月30日 星期一
文字特效
複習完 Matrix 後突然想寫這個，花了 15 分鐘...
#include <windows.h> #include <stdio.h> #include <time.h> #include <stdlib.h> int main (/* daviddr 091129 */) { const int W = 80, H = 24, L = 12; const float PL = 2.5; float pos[W], step[W]; int x, y, c; char s[] = " "; srand (time(0)); SetConsoleTitle (L"Art"); HDC hdc = GetDC (FindWindow (0, L"Art")); SetBkColor (hdc,0); for (x=0; x<W; ++x) { pos [x] = rand() % int(H*PL); step[x] = float (rand() % 7 + 5) / 5; } for (;;Sleep(120)) for (x=0; x<W; ++x ) { for (y=0; y<H; ++y) { c = int (pos[x]) - y; c = ((0<c && c<L)? 255-c*255/L: 0); *s = 33 + rand() % 93; SetTextColor (hdc, c>225? 0xFFFFFF: RGB(0,c,0)); TextOutA (hdc, x*8, y*16, s, 1); } pos[x] += step[x]; if (pos[x] >= H*PL) { pos[x] = 0; step[x] = float (rand() % 7 + 5) / 5; } } }截圖：
2009年11月26日 星期四
stackoverflow Part 2
beautiful code 系列：
首先是廣為人知的 Carmack Q3:
float InvSqrt (float x) { float xhalf = 0.5f*x; int i = *(int*)&x; i = 0x5f3759df - (i>>1); x = *(float*)&i; x = x*(1.5f - xhalf*x*x); return x; }Haskell 著名示例：
qsort [] = [] qsort (x:xs) = qsort (filter (< x) xs) ++ [x] ++ qsort (filter (>= x) xs)談最佳化時常提及的 Duff's device，雖然只是 loop unrolling~
void strcpy (to, from, count) register char *to, *from; register int count; { int n = (count + 7) / 8; switch (count % 8) { case 0: do { *to = *from++; case 7: *to = *from++; case 6: *to = *from++; case 5: *to = *from++; case 4: *to = *from++; case 3: *to = *from++; case 2: *to = *from++; case 1: *to = *from++; } while (--n > 0); } }Evil C 裡的另一例：
switch (c&3) while ((c-=4)>=0) { foo(); case 3: foo(); case 2: foo(); case 1: foo(); case 0: }theycallhimtom post。給定 adjacency matrix，計算各個最短路徑：
for (int i = 0; i < adj.length; i++) for (int j = 0; j < adj.length; j++) for (int k = 0; k < adj.length; k++) adj[j][k] = min (adj[j][k], adj[j][i] + adj[i][k]);位元倒置：
unsigned int reverse (register unsigned int x) { x = ((x & 0xaaaaaaaa) >> 1) | ((x & 0x55555555) << 1); x = ((x & 0xcccccccc) >> 2) | ((x & 0x33333333) << 2); x = ((x & 0xf0f0f0f0) >> 4) | ((x & 0x0f0f0f0f) << 4); x = ((x & 0xff00ff00) >> 8) | ((x & 0x00ff00ff) << 8); return (x >> 16) | (x << 16); }ob-C Mandelbrot:
main(){float A,B,P,Q,X,Y,d;int i,D=80,n=3120;for(X=-2.,Y=-1.5,d=6./D;B=2*A
*B+Y,A=P-Q+X,n;((P=A*A)+(Q=B*B)>4||++i>D)&&putchar(*((n--%D?X+=d/2,i<D?i%11
:11:(X=-2.0,Y+=d,12))+"Mandelbrot! \n"))&&(A=B=P=Q=i=0));}
Another one#include <unistd.h>; float o=0.075,h=1.5,T,r,O,l,I;int _,L=80,s=3200;main(){for(;s%L|| (h-=o,T= -2),s;4 -(r=O*O)<(l=I*I)|++ _==L&&write(1,(--s%L?_<L?--_ %6:6:7)+"World! \n",1)&&(O=I=l=_=r=0,T+=o /2))O=I*2*O+h,I=l+T-r;}Ken Perlin 更短的..
main(k){float i,j,r,x,y=-16;while(puts(""),y++<15)for(x
=0;x++<84;putchar(" .:-;!/>)|&IH%*#"[k&15]))for(i=k=r=0;
j=r*r-i*i-2+x/25,i=2*r*i+y/10,j*j+i*i<11&&k++<111;r=j);}
2009年11月25日 星期三
stackoverflow.com
底下是改寫自 stackoverflow 的各種有趣寫法，
原始創意乃源出於不同作者：
Zeroing structs without memset:
FStruct s = {0};
或
FStruct s = {};
void search (char* query) { http://www.google.com/search?q=query return; }然後... Peter Mortensen 這樣寫到：
Fred* f = new(ram) Fred(); http://www.parashift.com/c++-faq-lite/dtors.html#faq-11.10 f->~Fred();同樣是 Peter Mortensen，身份變換 :-)
class A {}; struct B { A a; operator A&() { return a; } }; void func(A a) { } int main() { A a, c; B b; a = c; func (b); //yeah baby a = b; //gotta love this }Daniel Martin，Perl 運算式解析
sub e{($_)=@_;$n='( *-?[.\d]++ *)'; s:\($n\)|$n(.)$n:($1,0,$2*$4,$2+$4,0,$2-$4)[7&ord$3]//$2/$4:e?e($_):$_}gnibbler，Golfscript 文字模式 Sierpinski’s Triangle:
' /\ /__\ '4/{).+: ;.{ \ ++}%\{.+}%+~ ]}@~(*n*
惡作劇用：<td align="right" style="text-align: left">也是惡作劇用，fork bomb：
:(){ :|:& };:
放在 source code 最上方.. 然後執行它://&>/dev/null;x="${0%.*}";[ ! "$x" -ot "$0" ]||(rm -f "$x";cc -o "$x" "$0")&&"$x" $*;exit
In Quake 2 src:double normals[][] = { #include "normals.txt" };Multi-character constants:
enum state { stopped = 'STOP', running = 'RUN!', waiting = 'WAIT', };三元運算：
(a? b: c) = (b<c ? b: c);Patrick post，同樣是三元運算 (q = b ? x : y)
bool b; int x, y; int q = -!!b&x|-!b&y;算 PI:
double fldpi () { double pi; asm ("fldpi" : "=t" (pi)); return pi; }Static linked lists：
struct llist { int a; struct llist* next;}; #define cons(x,y) (struct llist[]) {{x,y}} struct llist *list = cons(1, cons(2, cons(3, cons(4, NULL))));整數求和:
#define sum(...) \ _sum(sizeof((int []){ __VA_ARGS__ })\ / sizeof(int), (int []){ __VA_ARGS__ }) int _sum (size_t count, int values[]) { int s = 0; while (count--) s += values[count]; return s; } int main () {printf("%i", sum (1,2,3));}上例的變形：
void f (int len, int *n, char c) {...} #define SIZEOF(a) sizeof(a)/sizeof*(a) #define NOPAREN(...) __VA_ARGS__ #define F1(arr,c) {int v[]={NOPAREN arr}; f(SIZEOF(v),v,c);} #define F2(c,...) {int v[]=__VA_ARGS__; f(SIZEOF(v),v,c);}再一種應用：
#include <algorithm>; #define MAX(x, ...) ({\ typeof(x) ra[] = {(x), __VA_ARGS__}; \ *std::max_element (&ra[0], &ra[sizeof ra/sizeof 0[ra]]); \ }) int main() { int i = 12; std::cout << MAX (i+1,1,3,6,i); }Johannes Schaub 寫的：
template <template From, template To> union union_cast { From from; To to; union_cast (From from): from(from) { } To getTo() const { return to; } };這是真的嗎 :-)
const MyClass& x = MyClass(); // temporary exists as long as x is in scopeRube 累計遞增:
def foo(n) lambda { |i| n += i } end測驗題，不用空白來定義結構。By Bojan Resnik：
typedef/**/struct{const/**/char/**/c;unsigned/**/u;}nospace;vividos 寫的：
int Func() try { // if success return 0; } catch(...) { return -1; }goto 和 dtor 的關係：
struct A { ~A() {cout<<"dtor";} }; int main() { { A a; goto L; } L: cout<<"After L"; }Evil C：
int logical_xor (int a, int b) { return !a != !b; }From Canadian Mind Products:
char *p; switch (n) { case 1: p = "one"; if (0) case 2: p = "two"; if (0) case 3: p = "three"; printf ("%s", p); break; }自動初始化 a 中元素為 0（那麼，本來會做嗎？）:
class C { int a[10]; C(): a() {} };template bitfields，by Kaz Dragon:
template <size_t X, size_t Y> struct bitfield { char left : X; char right : Y; };出自 Comptrol：
template <typename T> class Creator { friend void appear() { // a new function ::appear(), but it doesn't … // exist until Creator is instantiated } }; Creator<void> miracle; // ::appear() is created at this point Creator<double> oops; // ERROR: ::appear() is created a second time! template <typename T> class Creator { friend void feed(Creator<T>*){ // every T generates a different … // function ::feed() } }; Creator<void> one; // generates ::feed(Creator<void>*) Creator<double> two; // generates ::feed(Creator<double>*)用於 non-tmpl 的類 boost::identity。
出自 Johannes Schaub：
template<typename T> struct id { typedef T type; }; // void (*f)(); // same id<void()>::type *f; // void (*f(void(*p)()))(int); // same id<void(int)>::type *f(id<void()>::type *p); // int (*p)[2] = new int[10][2]; // same id<int[2]>::type *p = new int[10][2]; // void (C::*p)(int) = 0; // same id<void(int)>::type C::*p = 0;原載於 Sumant Tambe 的 Blog:
template<class T> // (a) a base template void f( T ) { std::cout << "f(T)\n"; } template<> void f<>(int*) { // (b) an explicit specialization std::cout << "f(int *) specilization\n"; } template<class T> // (c) another, overloads (a) void f( T* ) { std::cout << "f(T *)\n"; } template<> void f<>(int*) { // (d) another identical explicit specialization std::cout << "f(int *) another specilization\n"; }用複雜方法做簡單的事? 出自 Marcus Lindblom:
template<int> class foo; template<> struct foo<0> { int* get() { return array; } int* array; }; template<int i> struct foo: public foo<i-1> { int* get() { return array + 1; } };
2009年11月11日 星期三
2009年11月8日 星期日
資料統計圖
Windows 下有好幾種產生資料統計圖的方法：
(1) 用 Excel
#include <ole2.h> ..... CLSID clsid; IDispatch pApp; CoInitialize (0); CLSIDFromProgID (L"Excel.Application", &clsid); CoCreateInstance (clsid, 0, CLSCTX_LOCAL_SERVER, IID_IDispatch, (void**)&pApp);(2) 用 wgunplot
將 data 存入 data.txt，繪圖命令存入 script.txt，
再調用：system ("wgunplot script.txt");
(3) 用 plplot 直接畫
http://plplot.sourceforge.net/examples.php
(4) R
http://www.r-project.org/
專用於 Statistical Computing。
這幾種方案中，(2) 的彈性、易用性最高。
若需要更多的統計計算，可採用 R。
2009年11月1日 星期日
交換變數值
以下是一些交換變數值的方法。
C++ 中交換 2 變數常使用：std::swap
若要自己設計 swap 的話，
一種簡單、正確、且速度最快的寫法是：
C 語言中，若想支援各種類型的變數交換，
一個 (不是很好的) 方法是：
除此之外，就只能老實將型別寫出來：
其實多浪費了兩個堆疊空間 (用來傳參數)。
若僅針對 4Bytes 型別做交換 (如char、float 與 int)，
又想寫個「通用」的 swap，一種有趣的寫法是：
他的未最佳化組語有點長：
不調用函式的情況下，針對 1~4Bytes 型別，
非 MMX 之 x86 組語最快的寫法或許是：
(因為 xchg 指令週期一般比 mov 慢好幾倍)
理論上是...... 沒有。但存在一些限制性的做法：
當 a 和 b 指向相同位址時，以上均得出 0。
(2)~(5) 都有可能發生 over/under flow，
此外，(5) 還有除零導致程式掛掉的問題。
交換「整個」陣列時，一般會先用指標代表原先陣列，
因此，只需交換「代表它」的指標，就能達成陣列的交換。
即使陣列長度不同，也能成功「交換」。
例如：
於是在 C 中可以這樣設計：
於是交換 2 個 Widget 時，只要交換其中的物件指標，便可達成物件交換。
把單個或多個元素，從某個起始位址「整批」移到另個位址，
有 2 種情況：
陣列間元素的「部分」抽換，使用 memcpy。
陣列內元素的「自身抽換」，使用 memmove。
自身抽換之元素位於頭尾時，便形成自旋：
C++ 中交換 2 變數常使用：std::swap
若要自己設計 swap 的話，
一種簡單、正確、且速度最快的寫法是：
template <typename T> inline void swap (T& a, T& b) { T t(a); a=b; b=t; }
C 語言中，若想支援各種類型的變數交換，
一個 (不是很好的) 方法是：
#define SWAP(T,a,b) {T t=a; a=b; b=t;}
驗證：int main () { float a=-88, b=57; SWAP (float, a, b); printf ("\n a=%f b=%f", a, b); char *c="alcohol", *d="water"; SWAP (char*, c, d); printf ("\n c=%s d=%s", c, d); }
除此之外，就只能老實將型別寫出來：
void swap_int (int* a, int* b)
{
int c = *a; *a = *b; *b = c;
}
不過上面的寫法在 compiler 未最佳化前，其實多浪費了兩個堆疊空間 (用來傳參數)。
若僅針對 4Bytes 型別做交換 (如char、float 與 int)，
又想寫個「通用」的 swap，一種有趣的寫法是：
void swap2 (void* a, void* b) { int t = *(int*)a; *(int*)a = *(int*)b; *(int*)b = t; }
他的未最佳化組語有點長：
int t = *(int*)a; mov eax,dword ptr [a]
mov ecx,dword ptr [eax]
mov dword ptr [t],ecx
*(int*)a = *(int*)b; mov eax,dword ptr [a]
mov ecx,dword ptr [ b ]
mov edx,dword ptr [ecx]
mov dword ptr [eax],edx
*(int*)b = t; mov eax,dword ptr [ b ]
mov ecx,dword ptr [t]
mov dword ptr [eax],ecx
驗證：int main () { float a=-88.5, b=57.12; swap2 (&a, &b); printf ("\n a=%f b=%f", a, b); char c='5', d='9'; swap2 (&c, &d); printf ("\n c=%c d=%c", c, d); }
不調用函式的情況下，針對 1~4Bytes 型別，
非 MMX 之 x86 組語最快的寫法或許是：
__asm {
mov eax, a
mov ebx, b
mov b , eax
mov a , ebx
}比較短，但慢些的寫法：(因為 xchg 指令週期一般比 mov 慢好幾倍)
_asm {
mov eax, a
xchg eax, b
mov a , eax
}同樣短，但更慢：
_asm {
xchg a, eax
xchg b, eax
xchg a, eax
}比較慢，又額外浪費堆疊空間的寫法是：
_asm{
push a
push b
pop a
pop b
}或是浮點數專用的：
_asm{
fld a
fld b
fstp a
fstp b
}
是否存在不使用額外變數，而能將變數值正確交換的方法呢？理論上是...... 沒有。但存在一些限制性的做法：
a ^= b ^= a ^= b; (1)
a -= b = (a+=b) - b; (2)
a -= b += a -= b = -b; (3)
a = a+b, b = a-b, a = a-b; (4)
a *= b, b = a/b, a /= b; (5)
(2)(3)(4) 其實是相同的邏輯。 當 a 和 b 指向相同位址時，以上均得出 0。
(2)~(5) 都有可能發生 over/under flow，
此外，(5) 還有除零導致程式掛掉的問題。
交換「整個」陣列時，一般會先用指標代表原先陣列，
因此，只需交換「代表它」的指標，就能達成陣列的交換。
即使陣列長度不同，也能成功「交換」。
例如：
//原始陣列 int A[] = { ... }; int B[] = { ... }; //用指標代表某個陣列 int * a = A; int * b = B; //交換所代表的陣列 int *t = a; a = b; b = t;
於是在 C 中可以這樣設計：
void swap (void** a, void** b) { void* t = *a; *a = *b; *b = t; } void show (void* a, int n) { int i, *p = (int*) a; //此處強制定成 int for (i=0; i<n; ++i) printf (" %d", p[i]); puts (""); } int main() { const int NUM = 8; int i; int A[NUM] = {1,2,3,4,5,6,7,8}; int B[NUM] = {11,12,13,14,15,16,17,18}; void* a = A; void* b = B; show (a, NUM); show (b, NUM); swap (&a, &b); show (a, NUM); show (b, NUM); }更一般化地，我們可用 Widget 模式包覆在非基本類型上，
於是交換 2 個 Widget 時，只要交換其中的物件指標，便可達成物件交換。
把單個或多個元素，從某個起始位址「整批」移到另個位址，
有 2 種情況：
陣列間元素的「部分」抽換，使用 memcpy。
陣列內元素的「自身抽換」，使用 memmove。
自身抽換之元素位於頭尾時，便形成自旋：
void show (int* a, int n) { for (int i=0; i<n; ++i) printf (" %d", a[i]); puts (""); } void left_rotate (int* a, int n) { int t = a[0]; memmove (a, a+1, (n-1) * sizeof(int)); a[n-1] = t; } int main() { const int NUM = 8; int i=3, a[NUM] = {1,2,3,4,5,6,7,8}; while (i--) left_rotate (a, NUM), show (a, NUM); }
時間檢測
這是個 nano-second 精度的計時器，所做的僅是計算 MSR 差值。
在 multi-core processors 下，
各 CPU 的 TSC 同步率若不高，計時結果也會有誤差，
故建構 TimeProfile 時，便設法讓行程在單一 CPU 上跑。
設計上很簡便，但存在一些缺點：
1. 當他前幾次執行時，算出來的時間間隔會比較長。
因為 CPU 要花時間將算法的相關指令與資料 load 到 L2 甚至 L1 裡，
並調整諸如分支預測等機制的最佳模式。
2. Process 會因 context switch 被切出去，但 MSR 仍持續累加，
造成計時誤差。為避免此情形，可在檢測算法前，
將 process 與其中的 main thread 之 priority 提高，
待算完後再行回復：
DWORD oldProcess = GetPriorityClass (GetCurrentProcess());
DWORD oldThread = GetThreadPriority (GetCurrentThread());
SetPriorityClass (GetCurrentProcess(), REALTIME_PRIORITY_CLASS);
SetThreadPriority (GetCurrentThread(), THREAD_PRIORITY_TIME_CRITICAL);
...........[Algorithm]...........
SetThreadPriority (GetCurrentThread(), oldThread);
SetPriorityClass (GetCurrentProcess(), oldProcess);
但這步驟可能會讓其他 thread 無法順暢執行，
為保持 TimeProfile 的「輕鬆」簡便，便不加上這層考量。
3. 現今許多 CPU 為了省電節能，當處於閒置狀態會降低主頻，
如： Intel 的 EIST、AMD 的 C&Q，
故 TimeProfile 中額外定義了 set_freq，
以便在「需要」時重新讀取 CPU 頻率。
4. 需注意：真正的計時類別尚有許多細節待考慮，
本程式僅是個以簡便為主的粗糙設計。
#include <stdio.h> #include <windows.h> class TimeProfile { unsigned __int64 Hz; unsigned __int64 cycle() //讀取 64-bit 的 { //Model Specific Register (MSR) _asm { xor eax, eax _emit 0x0F //使用 cpuid 讓 CPU 不採無序執行指令 _emit 0xA2 //讓接下來的指令照流水線運行 _emit 0x0F //注意: 586 後的 IA32 才支援 RDTSC _emit 0x31 //結果置於 EDX:EAX } } public: void set_freq() { LARGE_INTEGER freq; QueryPerformanceFrequency (&freq); Hz = freq.QuadPart; } TimeProfile() { set_freq(); //若 multi-core 同步率高, 可不用加下面那行 SetThreadAffinityMask (GetCurrentThread(), 1); } double time; void start() {time = (double) cycle();} void end() {time = double (cycle()-time) /Hz;} };測試碼：
double cal_pi (int i = 50) { double pi = 2; for (; i>0; i--) pi = pi * double(i)/(2*i+1) + 2; return pi; } int main (/* daviddr 081211 */) { TimeProfile tp; double pi; for (int i=1; i<52; i+=5) { tp.start(); pi = cal_pi (i); tp.end(); printf ("\n[%2d] PI = %.15f Time: %.9f", i, pi, tp.time); } return system ("pause"); }從結果可以觀察到：第一次執行時，
因為 CPU 要做些「準備工作」，算出的時間便較長。
[ 1] PI = 2.666666666666667 Time: 0.000004494
[ 6] PI = 3.132156732156732 Time: 0.000001445
[11] PI = 3.141358472520136 Time: 0.000001507
[16] PI = 3.141586396037061 Time: 0.000001674
[21] PI = 3.141592479958224 Time: 0.000001756
[26] PI = 3.141592648659952 Time: 0.000001881
[31] PI = 3.141592653447636 Time: 0.000002006
[36] PI = 3.141592653585648 Time: 0.000002130
[41] PI = 3.141592653589671 Time: 0.000002255
[46] PI = 3.141592653589790 Time: 0.000002379
[51] PI = 3.141592653589793 Time: 0.000002504
2009年10月31日 星期六
陣列下界
為何c語言的陣列索引是由0開始?
我們很難揣測 Ritchie 當年的心意，我的想法是：
當時一些可自訂陣列上下界的語言，如 Fortran，
在計算陣列元素位址時，需先扣去下界。
由於當時組合語言尋址範圍有限，最差情況可能只有
+-32768，當下界太大 (>32768) 或太小 (<-32768) 時，
下界無法在 Compiler time 靜態議決使之轉成 0 以加快元素存取速度，
只能每次存取時，實實在在的減去下界，這會拖慢存取速度。
於是，有 2 種解法：
(1) 設下語法限制，限定上下界範圍。
(2) 直接規定下界為 0。
選擇 (2) 可令陣列與指標間，具備合適的語法對應：
由於 C 中，「無索引標記的陣列名稱」會自動轉成
「指向陣列第一個元素的指標常數」，如：
short a[1];
單用 a 時，其值為 &a[下界]，
亦為 &(*((a)+(下界))) 或 &(*((下界)+(a)))，
令下界 = 0，則如下形式為等「值」：
a == &a[0] == &(*(a+0)) == &(*(0+a))
其中，deref (*) 與 addrof (&) 可相消解:
&(*(a+0)) == &*(a+0) == a+0 == a
又因 a == &a[0]，消去 & 得：
*(a+0) == *a == a[0]
由於 *a 中的 a 視為「指向陣列第一個元素的指標常數」，
對其 dereference 必為「陣列第一個元素」，
而 *a 又相當於 *(a+0)，為滿足 *(a+x) ~ a[x]，
故令 x = 0。
當然，直接由 *(a+x) ~ a[x] 這個對應關係，
也可以推想出 x = 0 的合理性，但「含意」會比較少 :p
以「陣列指標」觀點來看 (pointer to array，
如：Type (*a)[n]) 則問題可簡化至指標。
由於我們允許，並限定： *(a+0) 等值於 *a
當 pointer 的 offset 起始值為 0 時，才能滿足此限定，
否則形成矛盾：*(a + low_boundary) == *a
訂閱：
文章 (Atom)
