From 22cdfdcca0212df64e3352a834bbacda5686d732 Mon Sep 17 00:00:00 2001
From: Samo Penic <samo.penic@fe.uni-lj.si>
Date: Thu, 24 Apr 2014 16:26:12 +0000
Subject: [PATCH] added 3 different calculations for Kc

---
 src/timestep.c  |   11 +++++++----
 src/tape        |    4 ++--
 src/shcomplex.c |   10 +++++-----
 src/shcomplex.h |    2 +-
 4 files changed, 15 insertions(+), 12 deletions(-)

diff --git a/src/shcomplex.c b/src/shcomplex.c
index 7afe27e..7d8c9d8 100644
--- a/src/shcomplex.c
+++ b/src/shcomplex.c
@@ -130,9 +130,9 @@
 }
 
 
-ts_double calculateKc(ts_vesicle *vesicle){
-    ts_int min=4;
-    ts_int max=vesicle->sphHarmonics->l-6;
+ts_double calculateKc(ts_vesicle *vesicle, ts_int lmin, ts_int lmax){
+    ts_int min=lmin;
+    ts_int max=lmax; //vesicle->sphHarmonics->l-3;
     ts_long i,j;
     ts_double retval, bval;
     gsl_matrix *A=gsl_matrix_alloc(max-min,2);
@@ -146,7 +146,7 @@
     for(i=min;i<max;i++){
             gsl_matrix_set(A, i-min,0,(ts_double)((i-1)*(i+2)));
             gsl_matrix_set(A, i-min,1,(ts_double)((i-1)*(i+2)*(i+1)*i));
-            fprintf(stderr,"%e %e\n", gsl_matrix_get(A,i-min,0), gsl_matrix_get(A,i-min,1));
+//            fprintf(stderr,"%e %e\n", gsl_matrix_get(A,i-min,0), gsl_matrix_get(A,i-min,1));
             bval=0.0;
             //average for m from 0..l (only positive m's)
             for(j=0;j<=i;j++){
@@ -155,7 +155,7 @@
                 bval=bval/(ts_double)vesicle->sphHarmonics->N/(ts_double)(i+1);
 
             gsl_vector_set(b,i-min,1.0/bval);
-            fprintf(stderr,"%e\n", 1.0/gsl_vector_get(b,i-min));
+//            fprintf(stderr,"%e\n", 1.0/gsl_vector_get(b,i-min));
     }
 //    fprintf(stderr,"b[2]=%e\n",gsl_vector_get(b,1));
     gsl_linalg_QR_decomp(A,tau);
diff --git a/src/shcomplex.h b/src/shcomplex.h
index 72f8030..fdd0eb0 100644
--- a/src/shcomplex.h
+++ b/src/shcomplex.h
@@ -5,5 +5,5 @@
 ts_spharm *complex_sph_init(ts_vertex_list *vlist, ts_uint l);
 ts_bool complex_sph_free(ts_spharm *sph);
 ts_bool calculateUlmComplex(ts_vesicle *vesicle);
-ts_double calculateKc(ts_vesicle *vesicle);
+ts_double calculateKc(ts_vesicle *vesicle, ts_int lmin, ts_int lmax);
 #endif
diff --git a/src/tape b/src/tape
index 27abf01..0801375 100644
--- a/src/tape
+++ b/src/tape
@@ -44,7 +44,7 @@
 
 ####### Program Control ############
 #how many MC sweeps between subsequent records of states to disk
-mcsweeps=250
+mcsweeps=2500
 #how many initial mcsweeps*inititer MC sweeps before recording to disk?
 inititer=0
 #how many records do you want on the disk iteration are there in a run?
@@ -52,7 +52,7 @@
 
 
 ###### Spherical harmonics ###########
-spherical_harmonics_coefficients=15
+spherical_harmonics_coefficients=21
 
 #shut up if we are using cluster!!!
 quiet=false
diff --git a/src/timestep.c b/src/timestep.c
index d5fd583..c0a8ef1 100644
--- a/src/timestep.c
+++ b/src/timestep.c
@@ -15,7 +15,7 @@
 
 ts_bool run_simulation(ts_vesicle *vesicle, ts_uint mcsweeps, ts_uint inititer, ts_uint iterations, ts_uint start_iteration){
 	ts_uint i, j,k;
-	ts_double r0,kc;
+	ts_double r0,kc1,kc2,kc3;
 	ts_double l1,l2,l3,volume=0.0,area=0.0,vmsr,bfsr, vmsrt, bfsrt;
 	ts_ulong epochtime;
 	FILE *fd1;
@@ -24,7 +24,7 @@
 	if(fd==NULL){
 		fatal("Cannot open statistics.csv file for writing",1);
 	}
-	fprintf(fd, "Epoch OuterLoop VertexMoveSucessRate BondFlipSuccessRate Volume Area lamdba1 lambda2 lambda3 Kc\n");
+	fprintf(fd, "Epoch OuterLoop VertexMoveSucessRate BondFlipSuccessRate Volume Area lamdba1 lambda2 lambda3 Kc(2-9) Kc(6-9) Kc(2-end)\n");
 	centermass(vesicle);
 	cell_occupation(vesicle);
 	if(start_iteration<inititer) ts_fprintf(stdout, "Starting simulation (first %d x %d MC sweeps will not be recorded on disk)\n", inititer, mcsweeps);
@@ -56,7 +56,10 @@
 			    calculateUlmComplex(vesicle);
 			    storeUlmComplex2(vesicle);
 			    saveAvgUlm2(vesicle);
-                kc=calculateKc(vesicle);
+                kc1=calculateKc(vesicle, 2,9);
+                kc2=calculateKc(vesicle, 6,9);
+                kc3=calculateKc(vesicle, 2,vesicle->sphHarmonics->l);
+            
 				fd1=fopen("state.dat","w");
 				fprintf(fd1,"%e %e\n",vesicle->volume, getR0(vesicle));
 				for(k=0;k<vesicle->vlist->n;k++){
@@ -71,7 +74,7 @@
 				fclose(fd1);
             }
 
-			fprintf(fd, "%lu %u %e %e %1.16e %1.16e %1.16e %1.16e %1.16e %1.16e\n",epochtime,i,vmsr,bfsr,volume, area,l1,l2,l3,kc);
+			fprintf(fd, "%lu %u %e %e %1.16e %1.16e %1.16e %1.16e %1.16e %1.16e %1.16e %1.16e\n",epochtime,i,vmsr,bfsr,volume, area,l1,l2,l3,kc1, kc2, kc3);
 		    fflush(fd);	
 		//	sprintf(filename,"timestep-%05d.pov",i-inititer);
 		//	write_pov_file(vesicle,filename);

--
Gitblit v1.9.3