X-Git-Url: https://hackdaworld.org/gitweb/?a=blobdiff_plain;f=nlsop.c;h=c4fc825f28886eb99c75f769dc8674681346b16f;hb=df812e3fbcc8bb37863672e5371ae40cdf28b43e;hp=7dc862ca67842929a3662f750e19dbb1c89c08da;hpb=e08a152c3adeed6ccad3b19d1a40a2685d700612;p=physik%2Fnlsop.git diff --git a/nlsop.c b/nlsop.c index 7dc862c..c4fc825 100644 --- a/nlsop.c +++ b/nlsop.c @@ -50,7 +50,10 @@ int usage(void) printf("-p \t pressure offset (default %f)\n",B_AP); printf("-A \t slope of linear c distribution (default %f)\n",A_CD); printf("-B \t linear c distribution offset (default %f)\n",B_CD); + /* printf("-C \t initial c concentration (default %d)\n",CC); + */ + printf("-D \t diffusion rate from cryst to amorph cells (default %f)\n",D_R); puts("-L \t load from file"); puts("-S \t save to file"); puts("-R \t read from random file"); @@ -78,31 +81,115 @@ int process_cell(d3_lattice *d3_l,u32 x,u32 y,u32 z,int r,double a,double b,int } } } - // p*=*conc; - // printf("debug: p * conc = %f\n",p); if(!(*thiz&AMORPH)) { if(get_rand(URAND_MAX)<=p) { MAKE_AMORPH(thiz); - printf("debug: c->a %d - %d\n",*t_c,*conc); - *t_c=*t_c+1-*conc; - } else *t_c+=1; + // *t_c=*t_c+1-*conc; + } // else *t_c+=1; } else { /* assume 1-p probability */ if(get_rand(URAND_MAX)>p) { MAKE_CRYST(thiz); - printf("debug: a->c %d - %d\n",*t_c,*conc); - *t_c=*t_c+1+*conc; - } else *t_c+=1; + // *t_c=*t_c+1+*conc; + } // else *t_c+=1; } + + *t_c+=1; + + return 1; +} + +int distrib_c(d3_lattice *d3_l,double d_r,double a,double b) +{ + u32 x,y,z; + int i,j,k,c; + int offset,off; + int carry; + + /* put one c ion somewhere in the lattice */ + x=get_rand(d3_l->max_x); + y=get_rand(d3_l->max_y); + z=get_rand_lgp(d3_l->max_z,a,b); + *(d3_l->extra+x+y*d3_l->max_x+z*d3_l->max_x*d3_l->max_y)+=1; + + /* diffusion in layer */ + for(i=0;imax_x;i++) + { + for(j=0;jmax_y;j++) + { + for(k=0;kmax_z;k++) + { + offset=i+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y; + /* case amorph */ + if(*(d3_l->status+offset)&AMORPH) + { + /* look at neighbours and move c ions */ + for(c=-1;c<=1;c++) + { + if(c!=0) + { + off=((i+d3_l->max_x+c)%d3_l->max_x)+j*d3_l->max_x+k*d3_l->max_x*d3_l->max_y; + /* case neighbour not amorph */ + if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(d_r*(*(d3_l->extra+off))); + /* case neighbour amorph */ + else carry=(*(d3_l->extra+off)-*(d3_l->extra+offset))/2; + if(carry!=0) + { + *(d3_l->extra+offset)+=carry; + *(d3_l->extra+off)-=carry; + } + } + } + for(c=-1;c<=1;c++) + { + if(c!=0) + { + off=i+((j+c+d3_l->max_y)%d3_l->max_y)*d3_l->max_x+k*d3_l->max_x*d3_l->max_y; + /* case neighbour not amorph */ + if(!(*(d3_l->status+off)&AMORPH)) carry=(int)(d_r*(*(d3_l->extra+off))); + /* case neighbour amorph */ + else carry=(*(d3_l->extra+off)-*(d3_l->extra+offset))/2; + if(carry!=0) + { + *(d3_l->extra+offset)+=carry; + *(d3_l->extra+off)-=carry; + } + } + } + } else + /* case not amorph */ + { + /* look at neighbours and move c ions */ + for(c=-1;c<=1;c++) + { + if(c!=0) + { + off=i+((j+c+d3_l->max_y)%d3_l->max_y)*d3_l->max_x+k*d3_l->max_x*d3_l->max_y; + /* case neighbour not amorph */ + if(!(*(d3_l->status+off)&AMORPH)) + { + carry=(*(d3_l->extra+off)-*(d3_l->extra+offset))/2; + if(carry!=0) + { + *(d3_l->extra+offset)+=carry; + *(d3_l->extra+off)-=carry; + } + } + } + } + } + } /* for z */ + } /* for y */ + } /* for x */ return 1; } -int distrib_c(d3_lattice *d3_l,int t_c,double a,double b) +int distrib_c_old(d3_lattice *d3_l,int t_c,double a,double b) { int i,j,k,total,area; double sum; @@ -126,17 +213,14 @@ int distrib_c(d3_lattice *d3_l,int t_c,double a,double b) } } temp=(int)((i+1)*a+b)*t_c/(sum*area_h[i]); - // if(temp) - // { - for(j=0;jstatus+(i*area)+j)&AMORPH)) { - if(!(*(d3_l->status+(i*area)+j)&AMORPH)) - { - *(d3_l->extra+(i*area)+j)=temp; - total+=temp; - } + *(d3_l->extra+(i*area)+j)=temp; + total+=temp; } - // } + } left=(int)(((i+1)*a+b)*t_c/sum)%area_h[i]; while(left) { @@ -167,7 +251,7 @@ int distrib_c(d3_lattice *d3_l,int t_c,double a,double b) return 1; } -int save_to_file(char *sf,d3_lattice *d3_l) +int save_to_file(char *sf,d3_lattice *d3_l,info *my_inf) { int sf_fd,c; @@ -181,6 +265,11 @@ int save_to_file(char *sf,d3_lattice *d3_l) puts("failed saving d3 lattice struct"); return -1; } + if(write(sf_fd,my_inf,sizeof(info))max_x*d3_l->max_y*d3_l->max_z; if(write(sf_fd,d3_l->status,c*sizeof(unsigned char))max_x*d3_l->max_y*d3_l->max_z; if((d3_l->status=(unsigned char*)malloc(c*sizeof(unsigned char)))==NULL) { @@ -241,8 +335,7 @@ int main(int argc,char **argv) { u32 max_x,max_y,max_z,x,y,z,x_c,y_c,z_c; int i,quit,escape,nowait; - double a_el,b_el,a_cd,b_cd,a_ap,b_ap; - int cc,steps,range,refresh; + int refresh; char s_file[MAX_CHARS]; char l_file[MAX_CHARS]; char r_file[MAX_CHARS]; @@ -253,22 +346,34 @@ int main(int argc,char **argv) char conc_txt[MAX_TXT]; char steps_txt[MAX_TXT]; char cc_txt[MAX_TXT]; + char a_txt[MAX_TXT]; + char s_txt[MAX_TXT]; + char ap_txt[MAX_TXT]; + char el_txt[MAX_TXT]; + char cd_txt[MAX_TXT]; + char r_txt[MAX_TXT]; + char ap2_txt[MAX_TXT]; + char cd2_txt[MAX_TXT]; + char el2_txt[MAX_TXT]; + char dr_txt[MAX_TXT]; char *arg_v[MAX_ARGV]; d3_lattice d3_l; + info my_info; max_x=X; max_y=Y; max_z=Z; - steps=STEPS; - range=RANGE; + my_info.steps=STEPS; + my_info.range=RANGE; refresh=REFRESH; - a_el=A_EL; - b_el=B_EL; - a_cd=A_CD; - b_cd=B_CD; - a_ap=A_AP; - b_ap=B_AP; - cc=CC; + my_info.a_el=A_EL; + my_info.b_el=B_EL; + my_info.a_cd=A_CD; + my_info.b_cd=B_CD; + my_info.a_ap=A_AP; + my_info.b_ap=B_AP; + my_info.cc=CC; + my_info.d_r=D_R; nowait=0; quit=0; escape=0; @@ -289,10 +394,10 @@ int main(int argc,char **argv) nowait=1; break; case 'a': - a_el=atof(argv[++i]); + my_info.a_el=atof(argv[++i]); break; case 'b': - b_el=atof(argv[++i]); + my_info.b_el=atof(argv[++i]); break; case 'x': max_x=atoi(argv[++i]); @@ -315,28 +420,33 @@ int main(int argc,char **argv) break; */ case 's': - steps=atoi(argv[++i]); + my_info.steps=atoi(argv[++i]); break; case 'd': refresh=atoi(argv[++i]); break; case 'r': - range=atoi(argv[++i]); + my_info.range=atoi(argv[++i]); break; case 'f': - a_ap=atof(argv[++i]); + my_info.a_ap=atof(argv[++i]); break; case 'p': - b_ap=atof(argv[++i]); + my_info.b_ap=atof(argv[++i]); break; case 'A': - a_cd=atof(argv[++i]); + my_info.a_cd=atof(argv[++i]); break; case 'B': - b_cd=atof(argv[++i]); + my_info.b_cd=atof(argv[++i]); break; + /* case 'C': - cc=atoi(argv[++i]); + my_info.cc=atoi(argv[++i]); + break; + */ + case 'D': + my_info.d_r=atof(argv[++i]); break; case 'L': strcpy(l_file,argv[++i]); @@ -358,6 +468,8 @@ int main(int argc,char **argv) y=max_y/2-1; z=max_z/2-1; + + if(!strcmp(r_file,"")) rand_init(NULL); else rand_init(r_file); @@ -379,22 +491,57 @@ int main(int argc,char **argv) memset(d3_l.extra,0,i*sizeof(int)); } else { - load_from_file(l_file,&d3_l); + load_from_file(l_file,&d3_l,&my_info); d3_lattice_init(&argc,argv,&d3_l,d3_l.max_x,d3_l.max_y,d3_l.max_z); } d3_event_init(&d3_l); + strcpy(a_txt,"args:"); + sprintf(s_txt,"steps: %d",my_info.steps); + sprintf(r_txt,"pressure range: %d",my_info.range); + sprintf(ap_txt,"pressure faktor: %.2f",my_info.a_ap); + sprintf(ap2_txt,"pressure offset: %.2f",my_info.b_ap); + sprintf(el_txt,"energy loss slope: %.2f",my_info.a_el); + sprintf(el2_txt,"energy loss offset: %.2f",my_info.b_el); + sprintf(cd_txt,"c distrib slope: %.2f",my_info.a_cd); + sprintf(cd2_txt,"c distrib offset: %.2f",my_info.b_cd); + sprintf(dr_txt,"diffusion rate: %.2f",my_info.d_r); + arg_v[1]=x_txt; + arg_v[2]=y_txt; + arg_v[3]=z_txt; + arg_v[4]=NULL; + arg_v[5]=status_txt; + arg_v[6]=conc_txt; + arg_v[7]=NULL; + arg_v[8]=NULL; + arg_v[9]=NULL; + arg_v[10]=steps_txt;; + arg_v[11]=cc_txt; + arg_v[12]=NULL; + arg_v[13]=NULL; + arg_v[14]=a_txt; + arg_v[15]=NULL; + arg_v[16]=s_txt; + arg_v[17]=r_txt; + arg_v[18]=ap_txt; + arg_v[19]=ap2_txt; + arg_v[20]=el_txt; + arg_v[21]=el2_txt; + arg_v[22]=cd_txt; + arg_v[23]=cd2_txt; + arg_v[24]=dr_txt; + if(!strcmp(l_file,"")) { i=0; - while((i