int moldyn_init(t_moldyn *moldyn,int argc,char **argv) {
- //int ret;
-
- //ret=moldyn_parse_argv(moldyn,argc,argv);
- //if(ret<0) return ret;
-
memset(moldyn,0,sizeof(t_moldyn));
rand_init(&(moldyn->random),NULL,1);
moldyn->vis.dim.z=z;
}
- printf("[moldyn] dimensions in A and A^2 respectively:\n");
+ printf("[moldyn] dimensions in A and A^3 respectively:\n");
printf(" x: %f\n",moldyn->dim.x);
printf(" y: %f\n",moldyn->dim.y);
printf(" z: %f\n",moldyn->dim.z);
return 0;
}
+/*
+ * creating lattice functions
+ */
+
int create_lattice(t_moldyn *moldyn,u8 type,double lc,int element,double mass,
u8 attr,u8 bnum,int a,int b,int c) {
return ret;
}
+/* fcc lattice init */
+int fcc_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) {
+
+ int count;
+ int i,j;
+ t_3dvec o,r,n;
+ t_3dvec basis[3];
+ double help[3];
+ double x,y,z;
+
+ x=a*lc;
+ y=b*lc;
+ z=c*lc;
+
+ if(origin) v3_copy(&o,origin);
+ else v3_zero(&o);
+
+ /* construct the basis */
+ for(i=0;i<3;i++) {
+ for(j=0;j<3;j++) {
+ if(i!=j) help[j]=0.5*lc;
+ else help[j]=.0;
+ }
+ v3_set(&basis[i],help);
+ }
+
+ v3_zero(&r);
+ count=0;
+
+ /* fill up the room */
+ r.x=o.x;
+ while(r.x<x) {
+ r.y=o.y;
+ while(r.y<y) {
+ r.z=o.z;
+ while(r.z<z) {
+ v3_copy(&(atom[count].r),&r);
+ atom[count].element=1;
+ count+=1;
+ for(i=0;i<3;i++) {
+ v3_add(&n,&r,&basis[i]);
+ if((n.x<x+o.x)&&
+ (n.y<y+o.y)&&
+ (n.z<z+o.z)) {
+ v3_copy(&(atom[count].r),&n);
+ count+=1;
+ }
+ }
+ r.z+=lc;
+ }
+ r.y+=lc;
+ }
+ r.x+=lc;
+ }
+
+ /* coordinate transformation */
+ help[0]=x/2.0;
+ help[1]=y/2.0;
+ help[2]=z/2.0;
+ v3_set(&n,help);
+ for(i=0;i<count;i++)
+ v3_sub(&(atom[i].r),&(atom[i].r),&n);
+
+ return count;
+}
+
+int diamond_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) {
+
+ int count;
+ t_3dvec o;
+
+ count=fcc_init(a,b,c,lc,atom,origin);
+
+ o.x=0.25*lc;
+ o.y=0.25*lc;
+ o.z=0.25*lc;
+
+ if(origin) v3_add(&o,&o,origin);
+
+ count+=fcc_init(a,b,c,lc,&atom[count],&o);
+
+ return count;
+}
+
int add_atom(t_moldyn *moldyn,int element,double mass,u8 bnum,u8 attr,
t_3dvec *r,t_3dvec *v) {
int link_cell_update(t_moldyn *moldyn) {
int count,i,j,k;
- int nx,ny,nz;
+ int nx,ny;
t_atom *atom;
t_linkcell *lc;
+ double x,y,z;
atom=moldyn->atom;
lc=&(moldyn->lc);
nx=lc->nx;
ny=lc->ny;
- nz=lc->nz;
+
+ x=moldyn->dim.x/2;
+ y=moldyn->dim.y/2;
+ z=moldyn->dim.z/2;
for(i=0;i<lc->cells;i++)
- list_destroy_f(&(moldyn->lc.subcell[i]));
+ list_destroy_f(&(lc->subcell[i]));
for(count=0;count<moldyn->count;count++) {
- i=(atom[count].r.x+(moldyn->dim.x/2))/lc->x;
- j=(atom[count].r.y+(moldyn->dim.y/2))/lc->y;
- k=(atom[count].r.z+(moldyn->dim.z/2))/lc->z;
+ i=((atom[count].r.x+(moldyn->dim.x/2))/lc->x);
+ j=((atom[count].r.y+(moldyn->dim.y/2))/lc->y);
+ k=((atom[count].r.z+(moldyn->dim.z/2))/lc->z);
list_add_immediate_f(&(moldyn->lc.subcell[i+j*nx+k*nx*ny]),
&(atom[count]));
}
if(moldyn->pt_scale&(T_SCALE_BERENDSEN|T_SCALE_DIRECT))
scale_velocity(moldyn,FALSE);
if(moldyn->pt_scale&(P_SCALE_BERENDSEN|P_SCALE_DIRECT))
-{
-printf("going to do p scale ...\n");
scale_volume(moldyn);
-printf("done\n");
-}
/* check for log & visualization */
if(e) {
v3_add(&(atom[i].r),&(atom[i].r),&delta);
v3_scale(&delta,&(atom[i].f),0.5*tau_square/atom[i].mass);
v3_add(&(atom[i].r),&(atom[i].r),&delta);
-//if(i==5) printf("v: %f %f %f\n",atom[i].r.x,(atom[i].r.x+moldyn->dim.x/2)/moldyn->lc.x,2*atom[i].r.x/moldyn->dim.x);
check_per_bound(moldyn,&(atom[i].r));
-//if(i==5) printf("n: %f %f %f\n",atom[i].r.x,(atom[i].r.x+moldyn->dim.x/2)/moldyn->lc.x,2*atom[i].r.x/moldyn->dim.x);
/* velocities */
v3_scale(&delta,&(atom[i].f),0.5*tau/atom[i].mass);
v3_add(&(atom[i].v),&(atom[i].v),&delta);
}
-//moldyn_bc_check(moldyn);
/* neighbour list update */
link_cell_update(moldyn);
* periodic boundayr checking
*/
-int check_per_bound(t_moldyn *moldyn,t_3dvec *a) {
+inline int check_per_bound(t_moldyn *moldyn,t_3dvec *a) {
double x,y,z;
t_3dvec *dim;
return 0;
}
-
-/*
- * lattice creation functions
- */
-
-/* fcc lattice init */
-int fcc_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) {
-
- int count;
- int i,j;
- t_3dvec o,r,n;
- t_3dvec basis[3];
- double help[3];
- double x,y,z;
-
- x=a*lc;
- y=b*lc;
- z=c*lc;
-
- if(origin) v3_copy(&o,origin);
- else v3_zero(&o);
-
- /* construct the basis */
- for(i=0;i<3;i++) {
- for(j=0;j<3;j++) {
- if(i!=j) help[j]=0.5*lc;
- else help[j]=.0;
- }
- v3_set(&basis[i],help);
- }
-
- v3_zero(&r);
- count=0;
-
- /* fill up the room */
- r.x=o.x;
- while(r.x<x) {
- r.y=o.y;
- while(r.y<y) {
- r.z=o.z;
- while(r.z<z) {
- v3_copy(&(atom[count].r),&r);
- atom[count].element=1;
- count+=1;
- for(i=0;i<3;i++) {
- v3_add(&n,&r,&basis[i]);
- if((n.x<x+o.x)&&
- (n.y<y+o.y)&&
- (n.z<z+o.z)) {
- v3_copy(&(atom[count].r),&n);
- count+=1;
- }
- }
- r.z+=lc;
- }
- r.y+=lc;
- }
- r.x+=lc;
- }
-
- /* coordinate transformation */
- help[0]=x/2.0;
- help[1]=y/2.0;
- help[2]=z/2.0;
- v3_set(&n,help);
- for(i=0;i<count;i++)
- v3_sub(&(atom[i].r),&(atom[i].r),&n);
-
- return count;
-}
-
-int diamond_init(int a,int b,int c,double lc,t_atom *atom,t_3dvec *origin) {
-
- int count;
- t_3dvec o;
-
- count=fcc_init(a,b,c,lc,atom,origin);
-
- o.x=0.25*lc;
- o.y=0.25*lc;
- o.z=0.25*lc;
-
- if(origin) v3_add(&o,&o,origin);
-
- count+=fcc_init(a,b,c,lc,&atom[count],&o);
-
- return count;
-}
projects / physik / posic.git / blobdiff