悲催的科学匠人 - 冷水's blog
天河-1A下编译安装 HDF5 和 CGNS
天河-1A下编译安装 HDF5 和 CGNS
天河的开发环境太单薄了,很多库都没有,只能自己编译,没法像ubuntu那样apt-get一下就ok了
在hdf5的页面下载 zlib szip,配置方式为
./configure -prefix /home/jack/lib
然后 make;make install
配置 hdf5
./configure --prefix=/vol-th/home/jack/lib --with-zlib=/vol-th/home/jack/lib --with-szlib=/vol-th/home/jack/lib --enable-fortran --enable-cxx
然后 make;make install
配置 CGNS
cgns的安装配置有点古怪,需要我们自己把zlib的库文件(*.a)和头文件(*.h)拷贝到/vol-th/home/jack/lib下
./configure --prefix=/vol-th/home/jack/lib --enable-64bit --with-hdf5=/vol-th/home/jack/lib --with-zlib=/vol-th/home/jack/lib --with-szip=/vol-th/home/jack/lib/szip-2.1/szip/lib/libsz.a
这里 --enable-64bit 得看实际需要来确定是否加上
然后需要修改make.defs中的两个设置。原始文件没有-ldl,这里必须自己加上
CLIBS = -lm -ldl
FLIBS = -ldl
然后 make;make install
ExStream最近的进展
实现了若干新的扩展模块
- Gamma-Reynolds转捩模型
- 任意周期对接边界
- 运行时参数控制,不暂停计算的情况下改变某些重要参数
- 基于一阶迎风格式的GMRES+ILU线性方程组求解,以及二阶格式的Jacobian组装
以前读书的时候,没有好好学习的几个领域
以前读书的时候,没有好好学习的几个领域:
- 数字信号处理与随机过程
- 线性系统与自动控制
- 流动稳定性分析
最近花了N多时间补课。光知道数值分析,不理解物理原理,其实也干不了更多事。
ExChimera搁了近三个月,好多东西搞忘了。查看svn的提交log,都不知道当初在干什么了。shit
kmean算法实现
具体描述见machine learning in action
main函数是一个二维点cluster测试
#include <math.h> #include <stdlib.h> #include <stdio.h> struct KMEAN{ int nvec,dim,ngroup; double *data, /* [nvec][dim] */ *mu, /* [ngroup][dim] */ *scale, *min; /* [dim] */ int *groups, /* nvec */ *count; /* ngroup] */ }; struct KMEAN* KMEAN_init(int nvec, int dim, int ngroup) { struct KMEAN* kmean = (struct KMEAN*) calloc(1,sizeof(struct KMEAN)); kmean->nvec = nvec; kmean->dim = dim; kmean->ngroup = ngroup; kmean->data = (double*) calloc(nvec*dim, sizeof(double)); kmean->mu = (double*) calloc(ngroup*dim, sizeof(double)); kmean->scale = (double*) calloc(dim, sizeof(double)); kmean->min = (double*) calloc(dim, sizeof(double)); //kmean->max = (double*) calloc(dim, sizeof(double)); kmean->groups = (int*) calloc(nvec, sizeof(int)); //kmean->mark = (int*) calloc(nvec, sizeof(int)); kmean->count = (int*) calloc(ngroup, sizeof(int)); } void KMEAN_free(struct KMEAN* kmean) { free(kmean->data); free(kmean->mu); free(kmean->min); //free(kmean->max); free(kmean->scale); free(kmean->groups); //free(kmean->mark); free(kmean->count); } inline static double dataij(struct KMEAN* kmean, int i, int j) { return kmean->data[i*kmean->dim + j]; } inline static double muij(struct KMEAN *kmean, int i, int j) { return kmean->mu[i*kmean->dim + j]; } inline static void setdataij(struct KMEAN* kmean, int i, int j, double v) { kmean->data[i*kmean->dim + j] = v; } inline static void setmuij(struct KMEAN *kmean, int i, int j, double v) { kmean->mu[i*kmean->dim + j] = v; } inline static void addmuij(struct KMEAN *kmean, int i, int j, double v) { kmean->mu[i*kmean->dim + j] += v; } void KMEAN_OutputPlt(struct KMEAN* kmean, const char* fname) { FILE *fp = fopen(fname,"w"); int n,d; if(fp==NULL) {puts("Can not open file for outputing tecplot"); return;} fprintf(fp,"ZONE T=\"DATA-%s\"\n", fname); for(n=0;n<kmean->nvec;n++) { for(d=0;d<kmean->dim;d++) fprintf(fp,"%11.4E ", dataij(kmean,n,d)); fprintf(fp,"%d\n", kmean->groups[n] ); } fprintf(fp,"ZONE T=\"MU-%s\"\n", fname); for(n=0;n<kmean->ngroup;n++) { for(d=0;d<kmean->dim;d++) fprintf(fp,"%11.4E ", muij(kmean,n,d)); fprintf(fp,"%d\n", n ); } fclose(fp); } void KMEAN_CalcScale(struct KMEAN *kmean) { int d,n; double dmax,dmin; // compute range of each dim for(d=0;d<kmean->dim;d++){ dmax = dmin = dataij(kmean,0,d); for(n=1;n<kmean->nvec;n++){ double x; x = dataij(kmean,n,d); if( dmax<x ) dmax = x; if( dmin>x ) dmin = x; } kmean->scale[d] = dmax - dmin; kmean->min[d] = dmin; //kmean->max[d] = dmax; //printf("Scaling: DIM %4d min=%11.4E max=%11.4E\n",d,dmin,dmax); } // randomly init mu for(n=0;n<kmean->ngroup;n++){ //printf("init mu Group %3d ",n); for(d=0;d<kmean->dim;d++){ setmuij(kmean, n,d, kmean->min[d] + 0.8*(rand()/(double)(RAND_MAX)) * kmean->scale[d] ); //printf("x_%3d=%11.4E ",d,muij(kmean,n,d)); } //puts(""); } for(n=0;n<kmean->nvec;n++) kmean->groups[n] = 0; for(d=0;d<kmean->dim;d++) kmean->scale[d] = 1.0/ (kmean->scale[d] * kmean->scale[d]); } static int WhichGroup(struct KMEAN* kmean, int i) { int g; double mindist; int mingroup; for(g=0;g<kmean->ngroup;g++){ double dist=0; int d; for(d=0;d<kmean->dim;d++){ double xid,mgd; xid = dataij(kmean,i,d); mgd = muij(kmean,g,d); dist += (xid-mgd) * (xid-mgd) * kmean->scale[d]; } dist = sqrt(dist); if(g==0) { mindist = dist; mingroup = 0; }else if(mindist>dist){ mindist = dist; mingroup = g; } } return mingroup; } static double KMEAN_error(struct KMEAN* kmean) { double err=0.0; int n; for(n=0;n<kmean->nvec;n++) { int d,g; g = kmean->groups[n]; for(d=0;d<kmean->dim;d++) { double dd = dataij(kmean,n,d) - muij(kmean,g,d); err += dd*dd*kmean->scale[d]; } } return err; } static void UpdateMu(struct KMEAN* kmean) { int n,d; for(n=0;n<kmean->ngroup;n++) for(d=0;d<kmean->dim;d++) setmuij(kmean,n,d,0.0); for(n=0;n<kmean->nvec;n++) for(d=0;d<kmean->dim;d++){ addmuij(kmean, kmean->groups[n], d, dataij(kmean,n,d) ); } for(n=0;n<kmean->ngroup;n++) for(d=0;d<kmean->dim;d++) setmuij(kmean,n,d, muij(kmean,n,d)/kmean->count[n] ); } static void KMEAN_sweep(struct KMEAN* kmean, double *err, int *changed) { int n; *changed = 0; for(n=0;n<kmean->ngroup;n++) kmean->count[n] = 0; for(n=0;n<kmean->nvec;n++){ int g; g = WhichGroup(kmean,n); if(kmean->groups[n] != g) (*changed)++; kmean->groups[n] = g; kmean->count[ kmean->groups[n] ] ++; //printf("Vec %3d is Group %3d\n",n,g); } UpdateMu(kmean); *err = KMEAN_error(kmean); } double KMEAN_cluster(struct KMEAN *kmean) { int nchanged,it=0; double err; KMEAN_CalcScale(kmean); do{ KMEAN_sweep(kmean, &err, &nchanged); //printf("%6d changed %6d err=%11.4E\n",it++,nchanged, err); }while(nchanged>0); } void KMEAN_bisect(struct KMEAN *kmean, double *err) { int g, worstg; double err0,minderr; char fname[1024]; err0 = *err; minderr = err0; sprintf(fname,"skmean-g%d.plt",kmean->ngroup); KMEAN_OutputPlt(kmean, fname); // 找到最烂的分组 for(g=0;g<kmean->ngroup;g++) { struct KMEAN* subset = KMEAN_init(kmean->count[g], kmean->dim, 2); int i,si=0; double cerr; // fill group[g] into subset for(i=0;i<kmean->nvec;i++){ if(kmean->groups[i]==g){ int d; for(d=0;d<kmean->dim;d++) setdataij(subset,si,d, dataij(kmean,i,d)); //subset->mark[si] = kmean->mark[i]; si++; } } cerr = KMEAN_cluster(subset); KMEAN_free(subset); if(minderr > err0-cerr) {minderr = err0-cerr; worstg = g;} } g = worstg; // 最烂的分组 { struct KMEAN* subset = KMEAN_init(kmean->count[g], kmean->dim, 2); int i,si=0,d; double cerr; // fill group[g] into subset for(i=0;i<kmean->nvec;i++){ if(kmean->groups[i]==g){ for(d=0;d<kmean->dim;d++) setdataij(subset,si,d, dataij(kmean,i,d)); //subset->mark[si] = kmean->mark[i]; si++; } } // 更新mu for(d=0;d<kmean->dim;d++) { setmuij(kmean, g, d, muij(subset, 0, d) ); setmuij(kmean, kmean->ngroup, d, muij(subset, 1, d) ); } kmean->ngroup++; KMEAN_free(subset); // 重新对kmean做一次分组 *err = KMEAN_cluster(kmean); } return; } void KMEAN_BISECT(struct KMEAN* kmean) { int it=0,dest_ngroup = kmean->ngroup; double err; if(dest_ngroup<2) return; kmean->ngroup = 2; err = KMEAN_cluster(kmean); do{ KMEAN_bisect(kmean, &err); printf("Step %6d err=%11.4E\n",++it, err); }while(kmean->ngroup<dest_ngroup); } int main(void) { struct KMEAN *kmean = KMEAN_init(80, 2, 4); double err; { FILE* fp = fopen("testSet.txt","r"); int n=0; for(n=0;n<80;n++){ double x,y; fscanf(fp,"%lf%lf", &x,&y); setdataij(kmean,n,0, x); setdataij(kmean,n,1, y); } fclose(fp); } /* 下面这一行是采用普通kmean */ //err = KMEAN_cluster(kmean); printf("err=%11.4E\n", err); /* 下面这一行是采用二分kmean */ KMEAN_BISECT(kmean); KMEAN_OutputPlt(kmean, "skmean-g4.plt"); KMEAN_free(kmean); }