#include <omp.h>
#endif
-#ifdef PTHREADS
+#if defined PTHREADS || defined VISUAL_THREAD
#include <pthread.h>
#endif
#undef PSE_NAME
#undef PSE_COL
+#ifdef PTHREADS
+/* global mutexes */
+pthread_mutex_t *amutex;
+pthread_mutex_t emutex;
+#endif
+
/*
* the moldyn functions
*/
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;i<moldyn->count;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));
*/
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;
void *ptr;
t_atom *atom;
char name[16];
+#ifdef PTHREADS
+ pthread_mutex_t *mutex;
+#endif
new=a*b*c;
count=moldyn->count;
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;
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:
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 */
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 */
}
/* 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)
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;i<a;i++) {
r.y=o.y;
for(j=0;j<b;j++) {
r.z=o.z;
for(k=0;k<c;k++) {
+ switch(p_type) {
+ case PART_INSIDE_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;
+ 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;
}
/* 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)
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;i<a;i++) {
r.z=o.z;
for(k=0;k<c;k++) {
/* first atom */
+ switch(p_type) {
+ case PART_INSIDE_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;
- 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;
}
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;
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;
}
*
*/
-/* helper / special functions */
-
-#ifdef PTHREADS
-void *write_save_file(void *ptr) {
-
- int fd;
- char dir[128];
- t_moldyn *moldyn;
-
- moldyn=ptr;
-
- 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);
- if(fd<0) perror("[moldyn] save fd open");
- else {
- write(fd,moldyn,sizeof(t_moldyn));
- write(fd,moldyn->atom,moldyn->count*sizeof(t_atom));
- }
-
- close(fd);
-
- pthread_exit(NULL);
-
- return 0;
-}
-#endif
-
/* start the integration */
int moldyn_integrate(t_moldyn *moldyn) {
t_3dvec momentum;
t_moldyn_schedule *sched;
t_atom *atom;
-#ifndef PTHREADS
int fd;
char dir[128];
-#endif
double ds;
double energy_scale;
struct timeval t1,t2;
//double tp;
-#ifdef PTHREADS
+#ifdef VISUAL_THREAD
u8 first,change;
pthread_t io_thread;
int ret;
moldyn->debug=0;
/* zero & init moldyn copy */
-#ifdef PTHREADS
+#ifdef VISUAL_THREAD
memset(&md_copy,0,sizeof(t_moldyn));
atom_copy=malloc(moldyn->count*sizeof(t_atom));
if(atom_copy==NULL) {
}
#endif
+#ifdef PTHREADS
+ printf("##################\n");
+ printf("# USING PTHREADS #\n");
+ printf("##################\n");
+#endif
/* tell the world */
printf("[moldyn] integration start, go get a coffee ...\n");
}
if(s) {
if(!(moldyn->total_steps%s)) {
-#ifdef PTHREADS
+ 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);
+ if(fd<0) perror("[moldyn] save fd open");
+ else {
+ write(fd,moldyn,sizeof(t_moldyn));
+ write(fd,moldyn->atom,
+ moldyn->count*sizeof(t_atom));
+ }
+ close(fd);
+ }
+ }
+ 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);
md_copy.atom=atom_copy;
memcpy(atom_copy,moldyn->atom,moldyn->count*sizeof(t_atom));
change=0;
- ret=pthread_create(&io_thread,NULL,write_save_file,&md_copy);
+ ret=pthread_create(&io_thread,NULL,visual_atoms,&md_copy);
if(ret) {
- perror("[moldyn] create save file thread\n");
+ perror("[moldyn] create visual atoms thread\n");
return -1;
}
#else
- 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);
- if(fd<0) perror("[moldyn] save fd open");
- else {
- write(fd,moldyn,sizeof(t_moldyn));
- write(fd,moldyn->atom,
- moldyn->count*sizeof(t_atom));
- }
- close(fd);
-#endif
- }
- }
- if(a) {
- if(!(moldyn->total_steps%a)) {
visual_atoms(moldyn);
+#endif
}
}
/* 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,
jtom,
ktom,
bc_ik|bc_ij);
+
#ifdef STATIC_LISTS
}
#elif LOWMEM_LISTS
return 0;
}
+#ifdef VISUAL_THREAD
+void *visual_atoms(void *ptr) {
+#else
int visual_atoms(t_moldyn *moldyn) {
+#endif
int i;
char file[128+64];
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;
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 */
atom[i].ekin);
// bonds between atoms
+#ifndef VISUAL_THREAD
process_2b_bonds(moldyn,&vb,visual_bonds_process);
+#endif
// boundaries
if(dim.x) {
close(vb.fd);
+#ifdef VISUAL_THREAD
+ pthread_exit(NULL);
+
+}
+#else
+
return 0;
}
+#endif
/*
* fpu cntrol functions
projects / physik / posic.git / blobdiff