Merge branch 'leadoff'
[physik/posic.git] / search_bonds.c
index d3c1199..8291109 100644 (file)
 #include "moldyn.h"
 #include "potentials/albe.h"
 
+typedef struct s_data {
+       double pm,len;
+       u8 type;
+} t_data;
+
 int usage(char *prog) {
 
        printf("\nusage:\n");
@@ -28,25 +33,63 @@ int usage(char *prog) {
        return -1;
 }
 
+int process(t_moldyn *moldyn,t_atom *ai,t_atom *aj,void *ptr,u8 bc) {
+
+       t_3dvec dist;
+       double d;
+       t_data *data;
+
+       data=ptr;
+
+       if(aj<ai)
+               return 0;
+
+       switch(data->type) {
+               case 'a':
+                       if(ai->brand!=0)
+                               return 0;
+                       if(aj->brand!=0)
+                               return 0;
+                               
+                       break;
+               case 'b':
+                       if(ai->brand!=1)
+                               return 0;
+                       if(aj->brand!=1)
+                               return 0;
+                       break;
+               default:
+                       if(ai->brand==aj->brand)
+                               return 0;
+                       break;
+       }
+
+       v3_sub(&dist,&(ai->r),&(aj->r));
+       check_per_bound(moldyn,&dist);
+       d=v3_norm(&dist);
+
+       if((d<=data->len+data->pm)&(d>=data->len-data->pm))
+               printf(" # atoms %d/%d %d/%d - %f\n",
+                      ai->tag,ai->brand,
+                      aj->tag,aj->brand,d);
+
+       return 0;
+}
+
 int main(int argc,char **argv) {
 
        t_moldyn moldyn;
-       t_atom *itom,*jtom;
-       int i,j;
        int ret;
-       t_list n[27];
-       t_list *this;
-       t_linkcell *lc;
-       t_3dvec dist;
-       double d,bondlen,bondpm;
+       t_data data;
 
        if(argc!=5) {
                usage(argv[0]);
                return -1;
        }
 
-       bondlen=atof(argv[3]);
-       bondpm=atof(argv[4]);
+       data.type=argv[2][0];
+       data.len=atof(argv[3]);
+       data.pm=atof(argv[4]);
 
        memset(&moldyn,0,sizeof(t_moldyn));
 
@@ -58,67 +101,12 @@ int main(int argc,char **argv) {
        }
 
        /* link cell init */
-       moldyn.cutoff=bondlen+bondpm;
+       moldyn.cutoff=data.len+data.pm;
        link_cell_init(&moldyn,VERBOSE);
-       lc=&(moldyn.lc);
 
-       /* analyzing ... */
-       for(i=0;i<moldyn.count;i++) {
-               itom=&(moldyn.atom[i]);
-               link_cell_neighbour_index(&moldyn,
-                                         (itom->r.x+moldyn.dim.x/2)/lc->x,
-                                         (itom->r.y+moldyn.dim.y/2)/lc->y,
-                                         (itom->r.z+moldyn.dim.z/2)/lc->z,
-                                         n);
-               for(j=0;j<27;j++) {
-                       this=&(n[j]);
-                       list_reset_f(this);
-
-                       if(this->start==NULL)
-                               continue;
-
-                       do {
-
-                               jtom=this->current->data;
-
-                               if(jtom<itom)
-                                       continue;
-
-                               switch(argv[2][0]) {
-                                       case 'a':
-                                               if(itom->brand!=0)
-                                                       continue;
-                                               if(jtom->brand!=0)
-                                                       continue;
-                                               break;
-                                       case 'b':
-                                               if(itom->brand!=1)
-                                                       continue;
-                                               if(jtom->brand!=1)
-                                                       continue;
-                                               break;
-                                       default:
-                                               if(itom->brand==jtom->brand)
-                                                       continue;
-                                               break;
-                               }
-
-                               v3_sub(&dist,&(itom->r),&(jtom->r));
-                               check_per_bound(&moldyn,&dist);
-                               d=v3_norm(&dist);
-
-                               if((d<=bondlen+bondpm)&(d>=bondlen-bondpm)) {
-
-                                       printf(" # atoms %d/%d %d/%d - %f\n",
-                                              itom->tag,itom->brand,
-                                              jtom->tag,jtom->brand,d);
-
-                               }
-
-                       } while(list_next_f(this)!=L_NO_NEXT_ELEMENT);
-               }
-       }
+       process_2b_bonds(&moldyn,&data,process);
 
+       /* analyzing ... */
        link_cell_shutdown(&moldyn);
 
        moldyn_free_save_file(&moldyn);