From: hackbard <hackbard>
Date: Fri, 15 Dec 2006 15:49:06 +0000 (+0000)
Subject: pot hopefully ok, doing timestep tests for 450C now ...
X-Git-Url: https://hackdaworld.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c0a8b254109929fba10795e644187c51742108a8;p=physik%2Fposic.git

pot hopefully ok, doing timestep tests for 450C now ...
---

diff --git a/moldyn.c b/moldyn.c
index f8bfd3f..5aa6d87 100644
--- a/moldyn.c
+++ b/moldyn.c
@@ -767,23 +767,23 @@ int moldyn_add_schedule(t_moldyn *moldyn,int runs,double tau) {
 	t_moldyn_schedule *schedule;
 
 	schedule=&(moldyn->schedule);
-	count=++(schedule->content_count);
+	count=++(schedule->total_sched);
 
-	ptr=realloc(moldyn->schedule.runs,count*sizeof(int));
+	ptr=realloc(schedule->runs,count*sizeof(int));
 	if(!ptr) {
 		perror("[moldyn] realloc (runs)");
 		return -1;
 	}
-	moldyn->schedule.runs=ptr;
-	moldyn->schedule.runs[count-1]=runs;
+	schedule->runs=ptr;
+	schedule->runs[count-1]=runs;
 
 	ptr=realloc(schedule->tau,count*sizeof(double));
 	if(!ptr) {
 		perror("[moldyn] realloc (tau)");
 		return -1;
 	}
-	moldyn->schedule.tau=ptr;
-	moldyn->schedule.tau[count-1]=tau;
+	schedule->tau=ptr;
+	schedule->tau[count-1]=tau;
 
 	return 0;
 }
@@ -806,16 +806,16 @@ int moldyn_set_schedule_hook(t_moldyn *moldyn,void *hook,void *hook_params) {
 
 int moldyn_integrate(t_moldyn *moldyn) {
 
-	int i,sched;
+	int i;
 	unsigned int e,m,s,v;
 	t_3dvec p;
-	t_moldyn_schedule *schedule;
+	t_moldyn_schedule *sched;
 	t_atom *atom;
 	int fd;
 	char dir[128];
 	double ds;
 
-	schedule=&(moldyn->schedule);
+	sched=&(moldyn->schedule);
 	atom=moldyn->atom;
 
 	/* initialize linked cell method */
@@ -852,12 +852,12 @@ int moldyn_integrate(t_moldyn *moldyn) {
 	moldyn->debug=0;
 
 	/* executing the schedule */
-	for(sched=0;sched<moldyn->schedule.content_count;sched++) {
+	for(sched->count=0;sched->count<sched->total_sched;sched->count++) {
 
 		/* setting amount of runs and finite time step size */
-		moldyn->tau=schedule->tau[sched];
+		moldyn->tau=sched->tau[sched->count];
 		moldyn->tau_square=moldyn->tau*moldyn->tau;
-		moldyn->time_steps=schedule->runs[sched];
+		moldyn->time_steps=sched->runs[sched->count];
 
 	/* integration according to schedule */
 
@@ -907,7 +907,7 @@ int moldyn_integrate(t_moldyn *moldyn) {
 				visual_atoms(&(moldyn->vis),moldyn->time,
 				             moldyn->atom,moldyn->count);
 				printf("\rsched: %d, steps: %d, debug: %d",
-				       sched,i,moldyn->debug);
+				       sched->count,i,moldyn->debug);
 				fflush(stdout);
 			}
 		}
@@ -918,8 +918,8 @@ int moldyn_integrate(t_moldyn *moldyn) {
 	}
 
 		/* check for hooks */
-		if(schedule->hook)
-			schedule->hook(moldyn,schedule->hook_params);
+		if(sched->hook)
+			sched->hook(moldyn,sched->hook_params);
 
 		/* get a new info line */
 		printf("\n");
@@ -1108,6 +1108,9 @@ int potential_force_calc(t_moldyn *moldyn) {
 		}
 
 	}
+#ifdef DEBUG
+printf("\n\n");
+#endif
 
 	return 0;
 }
@@ -1397,6 +1400,14 @@ int tersoff_mult_2bp(t_moldyn *moldyn,t_atom *ai,t_atom *aj,u8 bc) {
 	/* add forces of 2bp (ij, ji) contribution
 	 * dVij = dVji and we sum up both: no 1/2) */
 	v3_add(&(ai->f),&(ai->f),&force);
+#ifdef DEBUG
+if(ai==&(moldyn->atom[0])) {
+	printf("dVij, dVji (2bp) contrib:\n");
+	printf("%f | %f\n",force.x,ai->f.x);
+	printf("%f | %f\n",force.y,ai->f.y);
+	printf("%f | %f\n",force.z,ai->f.z);
+}
+#endif
 
 	/* energy 2bp contribution (ij, ji) is 0.5 f_r f_c ... */
 	moldyn->energy+=(0.5*f_r*f_c);
@@ -1488,6 +1499,14 @@ int tersoff_mult_post_2bp(t_moldyn *moldyn,t_atom *ai,t_atom *aj,u8 bc) {
 
 	/* add force */
 	v3_add(&(ai->f),&(ai->f),&force);
+#ifdef DEBUG
+if(ai==&(moldyn->atom[0])) {
+	printf("dVij (3bp) contrib:\n");
+	printf("%f | %f\n",force.x,ai->f.x);
+	printf("%f | %f\n",force.y,ai->f.y);
+	printf("%f | %f\n",force.z,ai->f.z);
+}
+#endif
 
 	/* add energy of 3bp sum */
 	moldyn->energy+=(0.5*f_c*b*f_a);
@@ -1516,6 +1535,14 @@ int tersoff_mult_post_2bp(t_moldyn *moldyn,t_atom *ai,t_atom *aj,u8 bc) {
 
 	/* add force */
 	v3_add(&(ai->f),&(ai->f),&force);
+#ifdef DEBUG
+if(ai==&(moldyn->atom[0])) {
+	printf("dVji (3bp) contrib:\n");
+	printf("%f | %f\n",force.x,ai->f.x);
+	printf("%f | %f\n",force.y,ai->f.y);
+	printf("%f | %f\n",force.z,ai->f.z);
+}
+#endif
 
 	return 0;
 }
@@ -1777,6 +1804,14 @@ int tersoff_mult_3bp(t_moldyn *moldyn,t_atom *ai,t_atom *aj,t_atom *ak,u8 bc) {
 		v3_scale(&temp2,&temp2,tmp*B*exp(-mu*d_jk)*f_c_jk*0.5);
 		v3_add(&(ai->f),&(ai->f),&temp2); /* -1 skipped in f_a calc ^ */
 		                                  /* scaled with 0.5 ^ */
+#ifdef DEBUG
+if(ai==&(moldyn->atom[0])) {
+	printf("dVjk (3bp) contrib:\n");
+	printf("%f | %f\n",temp2.x,ai->f.x);
+	printf("%f | %f\n",temp2.y,ai->f.y);
+	printf("%f | %f\n",temp2.z,ai->f.z);
+}
+#endif
 
 	}
 
diff --git a/moldyn.h b/moldyn.h
index 9493f5f..38e4999 100644
--- a/moldyn.h
+++ b/moldyn.h
@@ -57,7 +57,8 @@ typedef struct s_linkcell {
 
 /* moldyn schedule structure */
 typedef struct s_moldyn_schedule {
-	int content_count;
+	int count;
+	int total_sched;
 	int *runs;
 	double *tau;
 	int (*hook)(void *moldyn,void *hook);
diff --git a/run b/run
index 3a15543..5a3a768 100755
--- a/run
+++ b/run
@@ -1,6 +1,6 @@
 mkdir -p saves video
 ./clean $1
-./sic $1
+./sic $@
 if [ "$?" == "0" ]; then
 	#./perms
 	if [ "$1" ] ; then
diff --git a/sic.c b/sic.c
index f4cdd7d..35931e0 100644
--- a/sic.c
+++ b/sic.c
@@ -14,9 +14,8 @@
 int main(int argc,char **argv) {
 
 	/* check argv */
-	if(argc!=2) {
-		printf("[sic] error: arg1 (vis/log/save location) ");
-		printf("must be given!\n");
+	if(argc!=3) {
+		printf("[sic] usage: %s <logdir> <temperatur>\n",argv[0]);
 		return -1;
 	}
 
@@ -32,7 +31,7 @@ int main(int argc,char **argv) {
 	double tau;
 
 	/* testing location & velocity vector */
-	//t_3dvec r,v;
+	t_3dvec r,v;
 
 	/* values */
 	tau=1.0e-15;	/* delta t = 1 fs */
@@ -119,29 +118,30 @@ int main(int argc,char **argv) {
 	moldyn_bc_check(&md);
 
 	/* testing configuration */
-	//r.x=2.45/2;	v.x=0;
+	//r.x=2.95/2;	v.x=0;
 	//r.y=0;		v.y=0;
 	//r.z=0;		v.z=0;
 	//add_atom(&md,SI,M_SI,0,
-	//           ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB,
+	//           ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP,//|ATOM_ATTR_HB,
 	//           &r,&v);
-	//r.x=-2.45/2;	v.x=0;
+	//r.x=-2.95/2;	v.x=0;
 	//r.y=0;		v.y=0;
 	//r.z=0;		v.z=0;
 	//add_atom(&md,SI,M_SI,0,
-	//           ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP|ATOM_ATTR_HB,
+	//           ATOM_ATTR_1BP|ATOM_ATTR_2BP|ATOM_ATTR_3BP,//|ATOM_ATTR_HB,
 	//           &r,&v);
 
 	/* setting a nearest neighbour distance for the moldyn checks */
 	set_nn_dist(&md,0.25*sqrt(3.0)*LC_SI); /* diamond ! */
 
 	/* set temperature */
-	printf("[sic] setting temperature\n");
+	printf("[sic] setting temperature -> %f\n",273+atof(argv[2]));
 	//set_temperature(&md,273.0+1410.0);
-	set_temperature(&md,273.0+450.0);
+	//set_temperature(&md,273.0+450.0);
 	//set_temperature(&md,273.0);
 	//set_temperature(&md,1.0);
 	//set_temperature(&md,0.0);
+	set_temperature(&md,atof(argv[2])+273.0);
 
 	/* set pressure */
 	printf("[sic] setting pressure\n");
@@ -159,13 +159,17 @@ int main(int argc,char **argv) {
 
 	/* create the simulation schedule */
 	printf("[sic] adding schedule\n");
-	moldyn_add_schedule(&md,20001,1.0);
+	moldyn_add_schedule(&md,20000,.1);
+	moldyn_add_schedule(&md,10000,.2);
+	moldyn_add_schedule(&md,6667,.3);
+	moldyn_add_schedule(&md,5000,.4);
+	moldyn_add_schedule(&md,4001,.5);
 
 	/* activate logging */
 	printf("[sic] activate logging\n");
 	moldyn_set_log_dir(&md,argv[1]);
-	moldyn_set_log(&md,LOG_TOTAL_ENERGY,200);
-	moldyn_set_log(&md,VISUAL_STEP,200);
+	moldyn_set_log(&md,LOG_TOTAL_ENERGY,100);
+	moldyn_set_log(&md,VISUAL_STEP,100);
 
 	/*
 	 * let's do the actual md algorithm now