From cef6cace886115944d4abc04f37a6741c437bce1 Mon Sep 17 00:00:00 2001
From: Samo Penic <samo.penic@gmail.com>
Date: Sat, 08 Mar 2014 11:50:45 +0000
Subject: [PATCH] debugged seg_fault. Found error in io.h. However, it is unknown why bond_no is not equal to neigh_no

---
 src/io.c |   78 ++++++++++++++++++++++++++++++++++----
 1 files changed, 69 insertions(+), 9 deletions(-)

diff --git a/src/io.c b/src/io.c
index 9e278ec..6de9654 100644
--- a/src/io.c
+++ b/src/io.c
@@ -9,6 +9,7 @@
 #include <sys/types.h>
 #include <dirent.h>
 #include "initial_distribution.h"
+#include "poly.h"
 
 ts_bool print_vertex_list(ts_vertex_list *vlist){
 	ts_uint i;
@@ -107,11 +108,11 @@
         fprintf(fh," %.17E\t%.17E\t%.17E\t%.17E\t%.17E\t%u\n",
         vlist->vtx[i]->xk,vlist->vtx[i]->c,vlist->vtx[i]->energy,
         vlist->vtx[i]->energy_h, vlist->vtx[i]->curvature, 0);
-        for(j=0;j<vlist->vtx[i]->neigh_no;j++){
+        for(j=0;j<vlist->vtx[i]->bond_no;j++){
             fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length_dual);
         }
             fprintf(fh,"\n");
-        for(j=0;j<vlist->vtx[i]->neigh_no;j++){
+        for(j=0;j<vlist->vtx[i]->bond_no;j++){
             fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length);
         }
             fprintf(fh,"\n");
@@ -206,7 +207,7 @@
 	ts_vertex_list *vlist=vesicle->vlist;
 	ts_bond_list *blist=vesicle->blist;
 	ts_vertex **vtx=vlist->vtx;
-    ts_uint i;
+    ts_uint i,j;
     	char filename[255];
 	FILE *fh;
 
@@ -217,29 +218,69 @@
 		return TS_FAIL;
 	}
 	/* Here comes header of the file */
+
+	//find number of extra vtxs and bonds of polymeres
+	ts_uint monono=0, polyno=0;
+	ts_bool poly=0;
+	if(vesicle->poly_list!=NULL){
+		if(vesicle->poly_list->poly[0]!=NULL){
+		polyno=vesicle->poly_list->n;
+		monono=vesicle->poly_list->poly[0]->vlist->n;
+		poly=1;
+		}
+	}
 	fprintf(fh, "<?xml version=\"1.0\"?>\n<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n <UnstructuredGrid>\n");
-    fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n, blist->n);
+    fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+monono*polyno, blist->n+monono*polyno);
     fprintf(fh,"<PointData Scalars=\"scalars\">\n<DataArray type=\"Int64\" Name=\"scalars\" format=\"ascii\">");
    	for(i=0;i<vlist->n;i++){
 		fprintf(fh,"%u ",vtx[i]->idx);
     }
+	//polymeres
+	if(poly){
+		for(i=0;i<vesicle->poly_list->n;i++){
+			for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
+				fprintf(fh,"%u ", vesicle->poly_list->poly[i]->vlist->vtx[j]->idx);
+			}
+		}
+	}
 
     fprintf(fh,"</DataArray>\n</PointData>\n<CellData>\n</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n");
 	for(i=0;i<vlist->n;i++){
 		fprintf(fh,"%e %e %e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
+	}
+	//polymeres
+	if(poly){
+		for(i=0;i<vesicle->poly_list->n;i++){
+			for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
+				fprintf(fh,"%e %e %e\n", vesicle->poly_list->poly[i]->vlist->vtx[j]->x,vesicle->poly_list->poly[i]->vlist->vtx[j]->y, vesicle->poly_list->poly[i]->vlist->vtx[j]->z );
+			}
+		}
 	}
 
     fprintf(fh,"</DataArray>\n</Points>\n<Cells>\n<DataArray type=\"Int64\" Name=\"connectivity\" format=\"ascii\">");
 	for(i=0;i<blist->n;i++){
 			fprintf(fh,"%u %u\n",blist->bond[i]->vtx1->idx,blist->bond[i]->vtx2->idx);
 	}
+	//polymeres
+	if(poly){
+		for(i=0;i<vesicle->poly_list->n;i++){
+			for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
+				fprintf(fh,"%u %u\n", vesicle->poly_list->poly[i]->blist->bond[j]->vtx1->idx,vesicle->poly_list->poly[i]->blist->bond[j]->vtx2->idx);
+			}
+	//grafted bonds
+		fprintf(fh,"%u %u\n", vesicle->poly_list->poly[i]->grafted_vtx->idx, vesicle->poly_list->poly[i]->vlist->vtx[0]->idx);
+		}
+
+	}
+	
+
     fprintf(fh,"</DataArray>\n<DataArray type=\"Int64\" Name=\"offsets\" format=\"ascii\">");
-    for (i=2;i<blist->n*2+1;i+=2){
+    for (i=2;i<(blist->n+monono*polyno)*2+1;i+=2){
     fprintf(fh,"%u ",i);
     }
     fprintf(fh,"\n");
     fprintf(fh,"</DataArray>\n<DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">\n");
-     for (i=0;i<blist->n;i++){
+     for (i=0;i<blist->n+monono*polyno;i++){
         fprintf(fh,"3 ");
     }
 
@@ -295,17 +336,24 @@
 
 
 ts_vesicle *parsetape(ts_uint *mcsweeps, ts_uint *inititer, ts_uint *iterations){
-    long int nshell=17,ncxmax=60, ncymax=60, nczmax=60;  // THIS IS DUE TO CONFUSE BUG!
+    long int nshell=17,ncxmax=60, ncymax=60, nczmax=60, npoly=10, nmono=20, pswitch=0;  // THIS IS DUE TO CONFUSE BUG!
     char *buf=malloc(255*sizeof(char));
     long int brezveze0=1;
     long int brezveze1=1;
     long int brezveze2=1;
-    ts_double xk0=25.0, dmax=1.67,stepsize=0.15;
+    ts_double xk0=25.0, dmax=1.67,stepsize=0.15,kspring=800.0,pressure=0.0;
 	long int iter=1000, init=1000, mcsw=1000;
+
+
     cfg_opt_t opts[] = {
         CFG_SIMPLE_INT("nshell", &nshell),
+        CFG_SIMPLE_INT("npoly", &npoly),
+        CFG_SIMPLE_INT("nmono", &nmono),
         CFG_SIMPLE_FLOAT("dmax", &dmax),
         CFG_SIMPLE_FLOAT("xk0",&xk0),
+	CFG_SIMPLE_INT("pswitch",&pswitch),
+	CFG_SIMPLE_FLOAT("pressure",&pressure),
+	CFG_SIMPLE_FLOAT("k_spring",&kspring),
         CFG_SIMPLE_FLOAT("stepsize",&stepsize),
         CFG_SIMPLE_INT("nxmax", &ncxmax),
         CFG_SIMPLE_INT("nymax", &ncymax),
@@ -335,18 +383,30 @@
 	*inititer=init;
 	*mcsweeps=mcsw;
 	vesicle=initial_distribution_dipyramid(nshell,ncxmax,ncymax,nczmax,stepsize);
+	vesicle->poly_list=init_poly_list(npoly,nmono, vesicle->vlist);
+	vesicle->spring_constant=kspring;
+	poly_assign_spring_const(vesicle);
+	
+
     vesicle->nshell=nshell;
     vesicle->dmax=dmax*dmax;
     vesicle->bending_rigidity=xk0;
+    vtx_set_global_values(vesicle); //copies xk0 to every vertex
+
+
     vesicle->stepsize=stepsize;
     vesicle->clist->ncmax[0]=ncxmax;
     vesicle->clist->ncmax[1]=ncymax;
     vesicle->clist->ncmax[2]=nczmax;
     vesicle->clist->max_occupancy=8;
-
+	vesicle->pressure=pressure/vesicle->bending_rigidity;	//all energy contributions need to be divided by bending_rigidity!
+    vesicle->pswitch=pswitch;
     cfg_free(cfg);
 	free(buf);
   //  fprintf(stderr,"NSHELL=%u\n",vesicle->nshell);
+
+			
+
     return vesicle;
 
 }

--
Gitblit v1.9.3