From 8e98c9cbbddf1c9c852874f956b0dbe2f9a4a922 Mon Sep 17 00:00:00 2001
From: hackbard <hackbard>
Date: Wed, 29 Nov 2006 13:23:18 +0000
Subject: [PATCH] test

---
 moldyn.c | 23 ++++++++++++++++-------
 moldyn.h | 10 +++-------
 sic.c    |  7 +++++--
 3 files changed, 24 insertions(+), 16 deletions(-)

diff --git a/moldyn.c b/moldyn.c
index 505f496..2f274c2 100644
--- a/moldyn.c
+++ b/moldyn.c
@@ -41,6 +41,7 @@ int moldyn_init(t_moldyn *moldyn,int argc,char **argv) {
 
 int moldyn_shutdown(t_moldyn *moldyn) {
 
+	moldyn_log_shutdown(moldyn);
 	link_cell_shutdown(moldyn);
 	moldyn_log_shutdown(moldyn);
 	rand_close(&(moldyn->random));
@@ -382,15 +383,19 @@ t_3dvec get_total_p(t_moldyn *moldyn) {
 	return p_total;
 }
 
-double estimate_time_step(t_moldyn *moldyn,double nn_dist,double t) {
+double estimate_time_step(t_moldyn *moldyn,double nn_dist) {
 
 	double tau;
 
-	tau=0.05*nn_dist/(sqrt(3.0*K_BOLTZMANN*t/moldyn->atom[0].mass));
-	tau*=1.0E-9;
-	if(tau<moldyn->tau)
-		printf("[moldyn] warning: time step  (%f > %.15f)\n",
-		       moldyn->tau,tau);
+	/* nn_dist is the nearest neighbour distance */
+
+	if(moldyn->t==5.0) {
+		printf("[moldyn] i do not estimate timesteps below %f K!\n",
+		       MOLDYN_CRITICAL_EST_TEMP);
+		return 23.42;
+	}
+
+	tau=(0.05*nn_dist*moldyn->atom[0].mass)/sqrt(3.0*K_BOLTZMANN*moldyn->t);
 
 	return tau;	
 }
@@ -744,6 +749,7 @@ int potential_force_calc(t_moldyn *moldyn) {
 	/* reset energy */
 	moldyn->energy=0.0;
 
+printf("DEBUG: count = %d\n",count);
 	for(i=0;i<count;i++) {
 	
 		/* reset force */
@@ -755,7 +761,8 @@ int potential_force_calc(t_moldyn *moldyn) {
 
 		/* 2 body pair potential/force */
 		if(atom[i].attr&(ATOM_ATTR_2BP|ATOM_ATTR_3BP)) {
-		
+	
+printf("DEBUG: processing atom %d\n",i);
 			link_cell_neighbour_index(moldyn,
 				(atom[i].r.x+moldyn->dim.x/2)/lc->x,
 				(atom[i].r.y+moldyn->dim.y/2)/lc->y,
@@ -765,6 +772,7 @@ int potential_force_calc(t_moldyn *moldyn) {
 			countn=lc->countn;
 			dnlc=lc->dnlc;
 
+printf("DEBUG: countn = %d - dnslc = %d\n",countn,dnlc);
 			for(j=0;j<countn;j++) {
 
 				this=&(neighbour[j]);
@@ -783,6 +791,7 @@ int potential_force_calc(t_moldyn *moldyn) {
 
 					if((btom->attr&ATOM_ATTR_2BP)&
 					   (atom[i].attr&ATOM_ATTR_2BP))
+printf("DEBUG: calling func2b\n");
 						moldyn->func2b(moldyn,
 						               &(atom[i]),
 						               btom,
diff --git a/moldyn.h b/moldyn.h
index 58302fb..947467c 100644
--- a/moldyn.h
+++ b/moldyn.h
@@ -118,12 +118,6 @@ typedef struct s_moldyn {
 	t_random random;	/* random interface */
 } t_moldyn;
 
-#define MOLDYN_LVSTAT_TOTAL_E		0x01
-#define MOLDYN_LVSTAT_TOTAL_M		0x02
-#define MOLDYN_LVSTAT_SAVE		0x04
-#define MOLDYN_LVSTAT_VISUAL		0x08
-#define MOLDYN_LVSTAT_INITIALIZED	0x10
-
 #define MOLDYN_STAT_PBX			0x08	/* periodic boudaries in x */
 #define MOLDYN_STAT_PBY			0x10	/* y */
 #define MOLDYN_STAT_PBZ			0x20	/* and z direction */
@@ -231,6 +225,8 @@ typedef struct s_tersoff_mult_params {
 #define MOLDYN_CUTOFF			1.0e-9
 #define MOLDYN_RUNS			1000000
 
+#define MOLDYN_CRITICAL_EST_TEMP	5.0
+
 #define MOLDYN_INTEGRATE_VERLET		0x00
 #define MOLDYN_INTEGRATE_DEFAULT 	MOLDYN_INTEGRATE_VERLET
 
@@ -307,7 +303,7 @@ double get_e_pot(t_moldyn *moldyn);
 double get_total_energy(t_moldyn *moldyn);
 t_3dvec get_total_p(t_moldyn *moldyn);
 
-double estimate_time_step(t_moldyn *moldyn,double nn_dist,double t);
+double estimate_time_step(t_moldyn *moldyn,double nn_dist);
 
 int link_cell_init(t_moldyn *moldyn);
 int link_cell_update(t_moldyn *moldyn);
diff --git a/sic.c b/sic.c
index 9c98857..3924375 100644
--- a/sic.c
+++ b/sic.c
@@ -30,6 +30,9 @@ int main(int argc,char **argv) {
 	/* misc variables, mainly to initialize stuff */
 	t_3dvec r,v;
 
+	/* temperature */
+	double t;
+
 	/* initialize moldyn */
 	printf("[sic] moldyn init\n");
 	moldyn_init(&md,argc,argv);
@@ -91,12 +94,12 @@ int main(int argc,char **argv) {
 
 	/* create the simulation schedule */
 	printf("[sic] adding schedule\n");
-	moldyn_add_schedule(&md,1000,1.0e-15);
+	moldyn_add_schedule(&md,10000,1.0e-12);
 
 	/* activate logging */
 	printf("[sic] activate logging\n");
 	moldyn_set_log(&md,LOG_TOTAL_ENERGY,"saves/test-energy",100);
-	moldyn_set_log(&md,VISUAL_STEP,"saves/test-visual",50);
+	moldyn_set_log(&md,VISUAL_STEP,"saves/test-visual",100);
 
 	/*
 	 * let's do the actual md algorithm now
-- 
2.20.1