X-Git-Url: https://hackdaworld.org/gitweb/?a=blobdiff_plain;f=potentials%2Falbe_fast.c;h=79162b1f31c62b2d95d01cded216e9a2ffc1584d;hb=HEAD;hp=603606e6b9ce081e108f64d7b33f544d009e35d9;hpb=815b8841eb76564b2f90421d4eae180253453106;p=physik%2Fposic.git diff --git a/potentials/albe_fast.c b/potentials/albe_fast.c index 603606e..79162b1 100644 --- a/potentials/albe_fast.c +++ b/potentials/albe_fast.c @@ -23,7 +23,7 @@ #ifdef PTHREADS #include -#define MAX_THREADS 4 +#define MAX_THREADS 2 #endif #include "../moldyn.h" @@ -39,12 +39,12 @@ extern pthread_mutex_t emutex; * virial calculation */ -#define albe_v_calc(a,f,d) a->virial.xx+=f->x*d->x; \ - a->virial.yy+=f->y*d->y; \ - a->virial.zz+=f->z*d->z; \ - a->virial.xy+=f->x*d->y; \ - a->virial.xz+=f->x*d->z; \ - a->virial.yz+=f->y*d->z +#define albe_v_calc(a,f,d) (a)->virial.xx+=(f)->x*(d)->x; \ + (a)->virial.yy+=(f)->y*(d)->y; \ + (a)->virial.zz+=(f)->z*(d)->z; \ + (a)->virial.xy+=(f)->x*(d)->y; \ + (a)->virial.xz+=(f)->x*(d)->z; \ + (a)->virial.yz+=(f)->y*(d)->z #ifndef PTHREADS @@ -452,7 +452,11 @@ int albe_potential_force_calc(t_moldyn *moldyn) { } /* force contribution for atom i */ +#ifdef MATTONI + scale=-0.5*(f_c*(df_r-b*df_a)); // - in albe formalism +#else scale=-0.5*(f_c*(df_r-b*df_a)+df_c*(f_r-b*f_a)); // - in albe formalism +#endif v3_scale(&force,&(dist_ij),scale); v3_add(&(ai->f),&(ai->f),&force); @@ -461,7 +465,8 @@ int albe_potential_force_calc(t_moldyn *moldyn) { v3_add(&(jtom->f),&(jtom->f),&force); /* virial */ - virial_calc(ai,&force,&(dist_ij)); + albe_v_calc(ai,&force,&(dist_ij)); + //virial_calc(ai,&force,&(dist_ij)); #ifdef DEBUG if(moldyn->time>DSTART&&moldyn->timetime>DSTART&&moldyn->timef),&(ai->f),&force); /* derivative wrt k */ +#ifdef MATTONI + v3_scale(&tmp,&dcosdrk,fcdg); + v3_scale(&force,&tmp,pre_dzeta); +#else v3_scale(&force,&dist_ik,-1.0*dfcg); // dri rik = - drk rik v3_scale(&tmp,&dcosdrk,fcdg); v3_add(&force,&force,&tmp); v3_scale(&force,&force,pre_dzeta); +#endif /* force contribution */ v3_add(&(ktom->f),&(ktom->f),&force); @@ -623,7 +634,8 @@ if(moldyn->time>DSTART&&moldyn->timedim.x/2)/lc->x, (itom[i].r.y+moldyn->dim.y/2)/lc->y, (itom[i].r.z+moldyn->dim.z/2)/lc->z, neighbour_i); - dnlc=lc->dnlc; - /* copy the neighbour lists */ #ifdef STATIC_LISTS #elif LOWMEM_LISTS @@ -827,7 +838,7 @@ void *potential_force_thread(void *ptr) { while(p!=-1) { jtom=&(itom[p]); - p=lc->subcell->list[p]; + p=lc->subcell->list[p]; // thread safe! #else this=&(neighbour_i[j]); list_reset_f(this); @@ -1080,22 +1091,33 @@ void *potential_force_thread(void *ptr) { } /* force contribution for atom i */ +#ifdef MATTONI + scale=-0.5*(f_c*(df_r-b*df_a)); // - in albe formalism +#else scale=-0.5*(f_c*(df_r-b*df_a)+df_c*(f_r-b*f_a)); // - in albe formalism +#endif v3_scale(&force,&(dist_ij),scale); - pthread_mutex_lock(&(amutex[ai->tag])); + if(pthread_mutex_lock(&(amutex[ai->tag]))) + perror("[albe fast] mutex lock (1)\n"); v3_add(&(ai->f),&(ai->f),&force); - pthread_mutex_unlock(&(amutex[ai->tag])); + if(pthread_mutex_unlock(&(amutex[ai->tag]))) + perror("[albe fast] mutex unlock (1)\n"); /* force contribution for atom j */ v3_scale(&force,&force,-1.0); // dri rij = - drj rij - pthread_mutex_lock(&(amutex[jtom->tag])); + if(pthread_mutex_lock(&(amutex[jtom->tag]))) + perror("[albe fast] mutex lock (2)\n"); v3_add(&(jtom->f),&(jtom->f),&force); - pthread_mutex_unlock(&(amutex[jtom->tag])); + if(pthread_mutex_unlock(&(amutex[jtom->tag]))) + perror("[albe fast] mutex unlock (2)\n"); /* virial */ - pthread_mutex_lock(&(amutex[ai->tag])); - virial_calc(ai,&force,&(dist_ij)); - pthread_mutex_unlock(&(amutex[ai->tag])); + if(pthread_mutex_lock(&(amutex[ai->tag]))) + perror("[albe fast] mutex lock (3)\n"); + albe_v_calc(ai,&force,&(dist_ij)); + //virial_calc(ai,&force,&(dist_ij)); + if(pthread_mutex_unlock(&(amutex[ai->tag]))) + perror("[albe fast] mutex unlock (3)\n"); #ifdef DEBUG if(moldyn->time>DSTART&&moldyn->timetime>DSTART&&moldyn->timeenergy+=energy; - pthread_mutex_unlock(&emutex); - pthread_mutex_lock(&(amutex[ai->tag])); + if(pthread_mutex_unlock(&emutex)) + perror("[albe fast] mutex unlock (energy)\n"); + if(pthread_mutex_lock(&(amutex[ai->tag]))) + perror("[albe fast] mutex lock (4)\n"); ai->e+=energy; - pthread_mutex_unlock(&(amutex[ai->tag])); + if(pthread_mutex_unlock(&(amutex[ai->tag]))) + perror("[albe fast] mutex unlock (4)\n"); /* reset k counter for second k loop */ kcount=0; @@ -1218,9 +1244,11 @@ if(moldyn->time>DSTART&&moldyn->timetag])); + if(pthread_mutex_lock(&(amutex[jtom->tag]))) + perror("[albe fast] mutex lock (5)\n"); v3_add(&(jtom->f),&(jtom->f),&force); - pthread_mutex_unlock(&(amutex[jtom->tag])); + if(pthread_mutex_unlock(&(amutex[jtom->tag]))) + perror("[albe fast] mutex unlock (5)\n"); #ifdef DEBUG if(moldyn->time>DSTART&&moldyn->timetime>DSTART&&moldyn->timetag])); - virial_calc(ai,&force,&dist_ij); + if(pthread_mutex_lock(&(amutex[ai->tag]))) + perror("[albe fast] mutex lock (6)\n"); + albe_v_calc(ai,&force,&dist_ij); + //virial_calc(ai,&force,&dist_ij); /* force contribution to atom i */ v3_scale(&force,&force,-1.0); v3_add(&(ai->f),&(ai->f),&force); - pthread_mutex_unlock(&(amutex[ai->tag])); + if(pthread_mutex_unlock(&(amutex[ai->tag]))) + perror("[albe fast] mutex unlock (6)\n"); /* derivative wrt k */ +#ifdef MATTONI + v3_scale(&tmp,&dcosdrk,fcdg); + v3_scale(&force,&tmp,pre_dzeta); +#else v3_scale(&force,&dist_ik,-1.0*dfcg); // dri rik = - drk rik v3_scale(&tmp,&dcosdrk,fcdg); v3_add(&force,&force,&tmp); v3_scale(&force,&force,pre_dzeta); +#endif /* force contribution */ - pthread_mutex_lock(&(amutex[ktom->tag])); + if(pthread_mutex_lock(&(amutex[ktom->tag]))) + perror("[albe fast] mutex lock (7)\n"); v3_add(&(ktom->f),&(ktom->f),&force); - pthread_mutex_unlock(&(amutex[ktom->tag])); + if(pthread_mutex_unlock(&(amutex[ktom->tag]))) + perror("[albe fast] mutex unlock (7)\n"); #ifdef DEBUG if(moldyn->time>DSTART&&moldyn->timetime>DSTART&&moldyn->timetag])); - virial_calc(ai,&force,&dist_ik); + if(pthread_mutex_lock(&(amutex[ai->tag]))) + perror("[albe fast] mutex lock (8)\n"); + albe_v_calc(ai,&force,&dist_ik); + //virial_calc(ai,&force,&dist_ik); /* force contribution to atom i */ v3_scale(&force,&force,-1.0); v3_add(&(ai->f),&(ai->f),&force); - pthread_mutex_unlock(&(amutex[ai->tag])); + if(pthread_mutex_unlock(&(amutex[ai->tag]))) + perror("[albe fast] mutex unlock (8)\n"); /* increase k counter */ kcount++; @@ -1368,7 +1409,7 @@ int albe_potential_force_calc(t_moldyn *moldyn) { /* prepare thread data */ pft_data[j].moldyn=moldyn; - pft_data[j].start=j*count/MAX_THREADS; + pft_data[j].start=j*(count/MAX_THREADS); if(j==MAX_THREADS-1) { pft_data[j].end=count; }