X-Git-Url: https://www.hackdaworld.org/gitweb/?a=blobdiff_plain;f=moldyn.c;h=33224a69f82fd62c175ea6e300f757f70b8851bb;hb=a29f0d8f1217a53ac4e1b90e65d24394f7c8967f;hp=a6f5c4557efc68b4ff01e6996792b782b5d2bcea;hpb=ff44bcba9df0d021568abcb553d1d490c442f696;p=physik%2Fposic.git

diff --git a/moldyn.c b/moldyn.c
index a6f5c45..33224a6 100644
--- a/moldyn.c
+++ b/moldyn.c
@@ -83,9 +83,6 @@ int moldyn_init(t_moldyn *moldyn,int argc,char **argv) {
 	pthread_mutex_init(&emutex,NULL);
 #endif
 
-	if(crtt)
-		printf("USING CRT\n");
-
 	return 0;
 }
 
@@ -2322,6 +2319,7 @@ int velocity_verlet(t_moldyn *moldyn) {
 
 		/* constraint relaxation */
 		if(crtt) {
+			// forces
 			basis_trafo(&(atom[i].f),FORWARD,
 			            trafo_angle[2*i],trafo_angle[2*i+1]);
 			if(constraints[3*i])
@@ -2332,6 +2330,17 @@ int velocity_verlet(t_moldyn *moldyn) {
 				atom[i].f.z=0;
 			basis_trafo(&(atom[i].f),BACKWARD,
 			            trafo_angle[2*i],trafo_angle[2*i+1]);
+			// velocities
+			basis_trafo(&(atom[i].v),FORWARD,
+			            trafo_angle[2*i],trafo_angle[2*i+1]);
+			if(constraints[3*i])
+				atom[i].v.x=0;
+			if(constraints[3*i+1])
+				atom[i].v.y=0;
+			if(constraints[3*i+2])
+				atom[i].v.z=0;
+			basis_trafo(&(atom[i].v),BACKWARD,
+			            trafo_angle[2*i],trafo_angle[2*i+1]);
 		}
 
 #ifndef QUENCH
@@ -2372,6 +2381,7 @@ int velocity_verlet(t_moldyn *moldyn) {
 
 		/* constraint relaxation */
 		if(crtt) {
+			// forces
 			basis_trafo(&(atom[i].f),FORWARD,
 			            trafo_angle[2*i],trafo_angle[2*i+1]);
 			if(constraints[3*i])
@@ -2382,6 +2392,17 @@ int velocity_verlet(t_moldyn *moldyn) {
 				atom[i].f.z=0;
 			basis_trafo(&(atom[i].f),BACKWARD,
 			            trafo_angle[2*i],trafo_angle[2*i+1]);
+			// velocities
+			basis_trafo(&(atom[i].v),FORWARD,
+			            trafo_angle[2*i],trafo_angle[2*i+1]);
+			if(constraints[3*i])
+				atom[i].v.x=0;
+			if(constraints[3*i+1])
+				atom[i].v.y=0;
+			if(constraints[3*i+2])
+				atom[i].v.z=0;
+			basis_trafo(&(atom[i].v),BACKWARD,
+			            trafo_angle[2*i],trafo_angle[2*i+1]);
 		}
 
 		/* again velocities [actually v(t+tau)] */