X-Git-Url: https://www.hackdaworld.org/gitweb/?a=blobdiff_plain;f=moldyn.c;h=8451c4436f43db2213bec83677c359fb4062f133;hb=787848dc9bf3daa2ccb76296de904688cc504a45;hp=3ac0a2d46b67c915bce1188b9707ead8aedd13b1;hpb=22a224e0f41f996954049c69bcc6fbb3cf62325f;p=physik%2Fposic.git diff --git a/moldyn.c b/moldyn.c index 3ac0a2d..8451c44 100644 --- a/moldyn.c +++ b/moldyn.c @@ -23,6 +23,10 @@ #include #endif +#if defined PTHREADS || defined VISUAL_THREAD +#include +#endif + #include "moldyn.h" #include "report/report.h" @@ -43,6 +47,12 @@ #undef PSE_NAME #undef PSE_COL +#ifdef PTHREADS +/* global mutexes */ +pthread_mutex_t *amutex; +pthread_mutex_t emutex; +#endif + /* * the moldyn functions */ @@ -62,13 +72,28 @@ int moldyn_init(t_moldyn *moldyn,int argc,char **argv) { rand_init(&(moldyn->random),NULL,1); moldyn->random.status|=RAND_STAT_VERBOSE; +#ifdef PTHREADS + pthread_mutex_init(&emutex,NULL); +#endif + return 0; } int moldyn_shutdown(t_moldyn *moldyn) { +#ifdef PTHREADS + int i; +#endif + printf("[moldyn] shutdown\n"); +#ifdef PTHREADS + for(i=0;icount;i++) + pthread_mutex_destroy(&(amutex[i])); + free(amutex); + pthread_mutex_destroy(&emutex); +#endif + moldyn_log_shutdown(moldyn); link_cell_shutdown(moldyn); rand_close(&(moldyn->random)); @@ -488,7 +513,8 @@ int moldyn_log_shutdown(t_moldyn *moldyn) { */ int create_lattice(t_moldyn *moldyn,u8 type,double lc,int element,double mass, - u8 attr,u8 brand,int a,int b,int c,t_3dvec *origin) { + u8 attr,u8 brand,int a,int b,int c,t_3dvec *origin, + u8 p_type,t_part_vals *p_vals) { int new,count; int ret; @@ -496,6 +522,9 @@ int create_lattice(t_moldyn *moldyn,u8 type,double lc,int element,double mass, void *ptr; t_atom *atom; char name[16]; +#ifdef PTHREADS + pthread_mutex_t *mutex; +#endif new=a*b*c; count=moldyn->count; @@ -518,6 +547,16 @@ int create_lattice(t_moldyn *moldyn,u8 type,double lc,int element,double mass, moldyn->atom=ptr; atom=&(moldyn->atom[count]); +#ifdef PTHREADS + ptr=realloc(amutex,(count+new)*sizeof(pthread_mutex_t)); + if(!ptr) { + perror("[moldyn] mutex realloc (add atom)"); + return -1; + } + amutex=ptr; + mutex=&(amutex[count]); +#endif + /* no atoms on the boundaries (only reason: it looks better!) */ if(!origin) { orig.x=0.5*lc; @@ -533,21 +572,21 @@ int create_lattice(t_moldyn *moldyn,u8 type,double lc,int element,double mass, switch(type) { case CUBIC: set_nn_dist(moldyn,lc); - ret=cubic_init(a,b,c,lc,atom,&orig); + ret=cubic_init(a,b,c,lc,atom,&orig,p_type,p_vals); strcpy(name,"cubic"); break; case FCC: if(!origin) v3_scale(&orig,&orig,0.5); set_nn_dist(moldyn,0.5*sqrt(2.0)*lc); - ret=fcc_init(a,b,c,lc,atom,&orig); + ret=fcc_init(a,b,c,lc,atom,&orig,p_type,p_vals); strcpy(name,"fcc"); break; case DIAMOND: if(!origin) v3_scale(&orig,&orig,0.25); set_nn_dist(moldyn,0.25*sqrt(3.0)*lc); - ret=diamond_init(a,b,c,lc,atom,&orig); + ret=diamond_init(a,b,c,lc,atom,&orig,p_type,p_vals); strcpy(name,"diamond"); break; default: @@ -575,6 +614,9 @@ int create_lattice(t_moldyn *moldyn,u8 type,double lc,int element,double mass, atom[ret].tag=count+ret; check_per_bound(moldyn,&(atom[ret].r)); atom[ret].r_0=atom[ret].r; +#ifdef PTHREADS + pthread_mutex_init(&(mutex[ret]),NULL); +#endif } /* update total system mass */ @@ -609,6 +651,16 @@ int add_atom(t_moldyn *moldyn,int element,double mass,u8 brand,u8 attr, moldyn->lc.subcell->list=ptr; #endif +#ifdef PTHREADS + ptr=realloc(amutex,(count+1)*sizeof(pthread_mutex_t)); + if(!ptr) { + perror("[moldyn] mutex realloc (add atom)"); + return -1; + } + amutex=ptr; + pthread_mutex_init(&(amutex[count]),NULL); +#endif + atom=moldyn->atom; /* initialize new atom */ @@ -659,12 +711,14 @@ int del_atom(t_moldyn *moldyn,int tag) { } /* cubic init */ -int cubic_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { +int cubic_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin, + u8 p_type,t_part_vals *p_vals) { int count; t_3dvec r; int i,j,k; t_3dvec o; + t_3dvec dist; count=0; if(origin) @@ -672,14 +726,39 @@ int cubic_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { else v3_zero(&o); + /* shift partition values */ + if(p_type) { + p_vals->p.x+=(a*lc)/2.0; + p_vals->p.y+=(b*lc)/2.0; + p_vals->p.z+=(c*lc)/2.0; + } + r.x=o.x; for(i=0;ip)); + if(v3_absolute_square(&dist)<(p_vals->r*p_vals->r)) { + v3_copy(&(atom[count].r),&r); + count+=1; + } + break; + case PART_OUTSIDE_R: + v3_sub(&dist,&r,&(p_vals->p)); + if(v3_absolute_square(&dist)>=(p_vals->r*p_vals->r)) { v3_copy(&(atom[count].r),&r); count+=1; + } + break; + default: + v3_copy(&(atom[count].r),&r); + count+=1; + break; + } r.z+=lc; } r.y+=lc; @@ -697,12 +776,14 @@ int cubic_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { } /* fcc lattice init */ -int fcc_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { +int fcc_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin, + u8 p_type,t_part_vals *p_vals) { int count; int i,j,k,l; t_3dvec o,r,n; t_3dvec basis[3]; + t_3dvec dist; count=0; if(origin) @@ -719,6 +800,13 @@ int fcc_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { basis[2].y=0.5*lc; basis[2].z=0.5*lc; + /* shift partition values */ + if(p_type) { + p_vals->p.x+=(a*lc)/2.0; + p_vals->p.y+=(b*lc)/2.0; + p_vals->p.z+=(c*lc)/2.0; + } + /* fill up the room */ r.x=o.x; for(i=0;ip)); + if(v3_absolute_square(&dist)<(p_vals->r*p_vals->r)) { v3_copy(&(atom[count].r),&r); count+=1; - r.z+=lc; + } + break; + case PART_OUTSIDE_R: + v3_sub(&dist,&r,&(p_vals->p)); + if(v3_absolute_square(&dist)>=(p_vals->r*p_vals->r)) { + v3_copy(&(atom[count].r),&r); + count+=1; + } + break; + default: + v3_copy(&(atom[count].r),&r); + count+=1; + break; + } /* the three face centered atoms */ for(l=0;l<3;l++) { v3_add(&n,&r,&basis[l]); + switch(p_type) { + case PART_INSIDE_R: + v3_sub(&dist,&n,&(p_vals->p)); + if(v3_absolute_square(&dist)<(p_vals->r*p_vals->r)) { + v3_copy(&(atom[count].r),&r); + count+=1; + } + break; + case PART_OUTSIDE_R: + v3_sub(&dist,&n,&(p_vals->p)); + if(v3_absolute_square(&dist)>=(p_vals->r*p_vals->r)) { + v3_copy(&(atom[count].r),&r); + count+=1; + } + break; + default: v3_copy(&(atom[count].r),&n); count+=1; + break; + } } + r.z+=lc; } r.y+=lc; } @@ -752,12 +876,13 @@ int fcc_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { return count; } -int diamond_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { +int diamond_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin, + u8 p_type,t_part_vals *p_vals) { int count; t_3dvec o; - count=fcc_init(a,b,c,lc,atom,origin); + count=fcc_init(a,b,c,lc,atom,origin,p_type,p_vals); o.x=0.25*lc; o.y=0.25*lc; @@ -765,7 +890,7 @@ int diamond_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) { if(origin) v3_add(&o,&o,origin); - count+=fcc_init(a,b,c,lc,&atom[count],&o); + count+=fcc_init(a,b,c,lc,&atom[count],&o,p_type,p_vals); return count; } @@ -1620,6 +1745,17 @@ int moldyn_integrate(t_moldyn *moldyn) { struct timeval t1,t2; //double tp; +#ifdef VISUAL_THREAD + u8 first,change; + pthread_t io_thread; + int ret; + t_moldyn md_copy; + t_atom *atom_copy; + + first=1; + change=0; +#endif + sched=&(moldyn->schedule); atom=moldyn->atom; @@ -1668,6 +1804,21 @@ int moldyn_integrate(t_moldyn *moldyn) { /* debugging, ignore */ moldyn->debug=0; + /* zero & init moldyn copy */ +#ifdef VISUAL_THREAD + memset(&md_copy,0,sizeof(t_moldyn)); + atom_copy=malloc(moldyn->count*sizeof(t_atom)); + if(atom_copy==NULL) { + perror("[moldyn] malloc atom copy (init)"); + return -1; + } +#endif + +#ifdef PTHREADS + printf("##################\n"); + printf("# USING PTHREADS #\n"); + printf("##################\n"); +#endif /* tell the world */ printf("[moldyn] integration start, go get a coffee ...\n"); @@ -1752,7 +1903,7 @@ int moldyn_integrate(t_moldyn *moldyn) { } if(s) { if(!(moldyn->total_steps%s)) { - snprintf(dir,128,"%s/s-%07.f.save", + snprintf(dir,128,"%s/s-%08.f.save", moldyn->vlsdir,moldyn->time); fd=open(dir,O_WRONLY|O_TRUNC|O_CREAT, S_IRUSR|S_IWUSR); @@ -1767,16 +1918,45 @@ int moldyn_integrate(t_moldyn *moldyn) { } if(a) { if(!(moldyn->total_steps%a)) { +#ifdef VISUAL_THREAD + /* check whether thread has not terminated yet */ + if(!first) { + ret=pthread_join(io_thread,NULL); + } + first=0; + /* prepare and start thread */ + if(moldyn->count!=md_copy.count) { + free(atom_copy); + change=1; + } + memcpy(&md_copy,moldyn,sizeof(t_moldyn)); + if(change) { + atom_copy=malloc(moldyn->count*sizeof(t_atom)); + if(atom_copy==NULL) { + perror("[moldyn] malloc atom copy (change)"); + return -1; + } + } + md_copy.atom=atom_copy; + memcpy(atom_copy,moldyn->atom,moldyn->count*sizeof(t_atom)); + change=0; + ret=pthread_create(&io_thread,NULL,visual_atoms,&md_copy); + if(ret) { + perror("[moldyn] create visual atoms thread\n"); + return -1; + } +#else visual_atoms(moldyn); +#endif } } /* display progress */ - if(!(moldyn->total_steps%10)) { + if(!(i%10)) { /* get current time */ gettimeofday(&t2,NULL); -printf("\rsched:%d, steps:%d/%d, T:%4.1f/%4.1f P:%4.1f/%4.1f V:%6.1f (%d)", +printf("sched:%d, steps:%d/%d, T:%4.1f/%4.1f P:%4.1f/%4.1f V:%6.1f (%d)\n", sched->count,i,moldyn->total_steps, moldyn->t,moldyn->t_avg, moldyn->p/BAR,moldyn->p_avg/BAR, @@ -2117,6 +2297,7 @@ int potential_force_calc(t_moldyn *moldyn) { jtom, ktom, bc_ik|bc_ij); + #ifdef STATIC_LISTS } #elif LOWMEM_LISTS @@ -2919,7 +3100,11 @@ int visual_bonds_process(t_moldyn *moldyn,t_atom *itom,t_atom *jtom, return 0; } +#ifdef VISUAL_THREAD +void *visual_atoms(void *ptr) { +#else int visual_atoms(t_moldyn *moldyn) { +#endif int i; char file[128+64]; @@ -2928,6 +3113,11 @@ int visual_atoms(t_moldyn *moldyn) { t_visual *v; t_atom *atom; t_vb vb; +#ifdef VISUAL_THREAD + t_moldyn *moldyn; + + moldyn=ptr; +#endif v=&(moldyn->vis); dim.x=v->dim.x; @@ -2937,17 +3127,19 @@ int visual_atoms(t_moldyn *moldyn) { help=(dim.x+dim.y); - sprintf(file,"%s/atomic_conf_%07.f.xyz",v->fb,moldyn->time); + sprintf(file,"%s/atomic_conf_%08.f.xyz",v->fb,moldyn->time); vb.fd=open(file,O_WRONLY|O_CREAT|O_TRUNC,S_IRUSR|S_IWUSR); if(vb.fd<0) { perror("open visual save file fd"); +#ifndef VISUAL_THREAD return -1; +#endif } /* write the actual data file */ // povray header - dprintf(vb.fd,"# [P] %d %07.f <%f,%f,%f>\n", + dprintf(vb.fd,"# [P] %d %08.f <%f,%f,%f>\n", moldyn->count,moldyn->time,help/40.0,help/40.0,-0.8*help); // atomic configuration @@ -2961,7 +3153,9 @@ int visual_atoms(t_moldyn *moldyn) { atom[i].ekin); // bonds between atoms +#ifndef VISUAL_THREAD process_2b_bonds(moldyn,&vb,visual_bonds_process); +#endif // boundaries if(dim.x) { @@ -3007,8 +3201,15 @@ int visual_atoms(t_moldyn *moldyn) { close(vb.fd); +#ifdef VISUAL_THREAD + pthread_exit(NULL); + +} +#else + return 0; } +#endif /* * fpu cntrol functions