From bbe7df469c40beb309542c01f25c262340978629 Mon Sep 17 00:00:00 2001
From: Samo Penic <samo.penic@fe.uni-lj.si>
Date: Thu, 15 Sep 2016 08:04:34 +0000
Subject: [PATCH] Merge branch 'pegs'

---
 src/main.c                 |    3 
 src/poly.c                 |  162 +++++++++++++++++++++++++++++++++++++++++------------
 src/io.c                   |    2 
 src/tape                   |    6 +-
 src/initial_distribution.c |    2 
 src/frame.c                |    5 +
 6 files changed, 135 insertions(+), 45 deletions(-)

diff --git a/src/frame.c b/src/frame.c
index d81ba7e..68b8a84 100644
--- a/src/frame.c
+++ b/src/frame.c
@@ -72,20 +72,23 @@
     }
 
 //Add all polymers to cells
+if(vesicle->poly_list!=NULL){
     for(i=0;i<vesicle->poly_list->n;i++){
 	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
     	cellidx=vertex_self_avoidance(vesicle, vesicle->poly_list->poly[i]->vlist->vtx[j]);
     	cell_add_vertex(vesicle->clist->cell[cellidx],vesicle->poly_list->poly[i]->vlist->vtx[j]);
 	}
     }
+}
 //Add all filaments to cells
+if(vesicle->filament_list!=NULL){
      for(i=0;i<vesicle->filament_list->n;i++){
 	for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){
     	cellidx=vertex_self_avoidance(vesicle, vesicle->filament_list->poly[i]->vlist->vtx[j]);
     	cell_add_vertex(vesicle->clist->cell[cellidx],vesicle->filament_list->poly[i]->vlist->vtx[j]);
 	}
     }
-   
+}   
 
     return TS_SUCCESS;
 }
diff --git a/src/initial_distribution.c b/src/initial_distribution.c
index 45a4124..71a2c4f 100644
--- a/src/initial_distribution.c
+++ b/src/initial_distribution.c
@@ -100,7 +100,7 @@
 	vesicle->clist->ncmax[0]=tape->ncxmax;
 	vesicle->clist->ncmax[1]=tape->ncymax;
 	vesicle->clist->ncmax[2]=tape->nczmax;
-	vesicle->clist->max_occupancy=8; /* hard coded max occupancy? */
+	vesicle->clist->max_occupancy=16; /* hard coded max occupancy? */
 
 	vesicle->pressure= tape->pressure;
 	vesicle->pswitch=tape->pswitch;
diff --git a/src/io.c b/src/io.c
index 2c76fd2..a478ca0 100644
--- a/src/io.c
+++ b/src/io.c
@@ -423,7 +423,7 @@
     vesicle->tape=parsetape(command_line_args.tape_fullfilename);
 // recreating space for cells // 
     vesicle->clist=init_cell_list(vesicle->tape->ncxmax, vesicle->tape->ncymax, vesicle->tape->nczmax, vesicle->tape->stepsize);
-	vesicle->clist->max_occupancy=8;
+	vesicle->clist->max_occupancy=16;
 //    vesicle->tape=(ts_tape *)malloc(sizeof(ts_tape));
 //    retval=fread(vesicle->tape, sizeof(ts_tape),1,fh);
 	retval=fread(iteration,sizeof(ts_uint),1,fh);
diff --git a/src/main.c b/src/main.c
index 596d4a5..cb6dd05 100644
--- a/src/main.c
+++ b/src/main.c
@@ -127,7 +127,8 @@
 		}
 	}
 			//printf("nucleus coords: %.17e %.17e %.17e\n",vesicle->nucleus_center[0], vesicle->nucleus_center[1], vesicle->nucleus_center[2]);
-
+//	write_vertex_xml_file(vesicle,0);
+//	exit(1);
 	run_simulation(vesicle, tape->mcsweeps, tape->inititer, tape->iterations, start_iteration);
 	write_master_xml_file(command_line_args.output_fullfilename);
 	write_dout_fcompat_file(vesicle,"dout");
diff --git a/src/poly.c b/src/poly.c
index 7c7450b..bcf01e5 100644
--- a/src/poly.c
+++ b/src/poly.c
@@ -6,7 +6,8 @@
 #include"bond.h"
 #include<math.h>
 #include"energy.h"
-
+#include"cell.h"
+#include"frame.h"
 ts_bool poly_assign_filament_xi(ts_vesicle *vesicle, ts_tape *tape){
 	ts_uint i;
 
@@ -47,45 +48,20 @@
 	poly->blist->bond[i]->bond_length=sqrt(vtx_distance_sq(poly->blist->bond[i]->vtx1,poly->blist->bond[i]->vtx2));
 	bond_energy(poly->blist->bond[i],poly);
 	}
-
+	vertex_list_assign_id(poly->vlist,TS_ID_FILAMENT);
 	return poly;
 }
 
 
-ts_poly_list *init_poly_list(ts_uint n_poly, ts_uint n_mono, ts_vertex_list *vlist, ts_vesicle *vesicle){
-	ts_poly_list *poly_list=(ts_poly_list *)calloc(1,sizeof(ts_poly_list));
-	poly_list->poly	= (ts_poly **)calloc(n_poly,sizeof(ts_poly *));
-	ts_uint i=0,j=0; //idx;
-	ts_uint gvtxi;
-	ts_double xnorm,ynorm,znorm,normlength;
-	ts_double dphi,dh;
 
-	// Grafting polymers:
-	if (vlist!=NULL){
-		if (n_poly > vlist->n){fatal("Number of polymers larger than numbero f vertices on a vesicle.",310);}
-	
-		while(i<n_poly){
-			gvtxi = rand() % vlist->n;
-			if (vlist->vtx[gvtxi]->grafted_poly == NULL){
-			poly_list->poly[i] = init_poly(n_mono, vlist->vtx[gvtxi]);
-			i++;
-			}
-		}
-	}
-	else
-	{
-		for(i=0;i<n_poly;i++){
-			poly_list->poly[i] = init_poly(n_mono, NULL);
-		}
-	}
-
-	poly_list->n = n_poly;
-
-	if (vlist!=NULL){
+ts_bool poly_initial_distribution(ts_poly_list *poly_list, ts_int i, ts_vesicle *vesicle){
 	/* Make straight grafted poylmers normal to membrane (polymer brush). Dist. between poly vertices put to 1*/
+		ts_double xnorm,ynorm,znorm,normlength;
 		ts_int intpoly=vesicle->tape->internal_poly;
-		for (i=0;i<poly_list->n;i++){
-	
+		ts_int cellidx;
+		ts_double posX,posY,posZ,prevPosX,prevPosY,prevPosZ, phi,costheta,sintheta;
+		ts_bool retval;
+		ts_int j,k,l,m;
 			xnorm=0.0;
 			ynorm=0.0;
 			znorm=0.0;
@@ -101,13 +77,123 @@
 			xnorm=xnorm/normlength;
 			ynorm=ynorm/normlength;
 			znorm=znorm/normlength;
-
+			//prepare starting position for building the polymeres
+			prevPosX=poly_list->poly[i]->grafted_vtx->x;
+			prevPosY=poly_list->poly[i]->grafted_vtx->y;
+			prevPosZ=poly_list->poly[i]->grafted_vtx->z;
 			for (j=0;j<poly_list->poly[i]->vlist->n;j++){
-				poly_list->poly[i]->vlist->vtx[j]->x = poly_list->poly[i]->grafted_vtx->x + xnorm*(ts_double)(j+1);
-				poly_list->poly[i]->vlist->vtx[j]->y = poly_list->poly[i]->grafted_vtx->y + ynorm*(ts_double)(j+1);
-				poly_list->poly[i]->vlist->vtx[j]->z = poly_list->poly[i]->grafted_vtx->z + znorm*(ts_double)(j+1);
+				//if(j==0){
+					posX=prevPosX+xnorm*(vesicle->clist->dmin_interspecies);
+					posY=prevPosY+ynorm*(vesicle->clist->dmin_interspecies);
+					posZ=prevPosZ+znorm*(vesicle->clist->dmin_interspecies);
+				//}else{
+				//	posX=prevPosX+xnorm;
+				//	posY=prevPosY+ynorm;
+				//	posZ=prevPosZ+znorm;
+				//}
+				//trying to go towards normal
+				k=0;
+				while(1){
+					poly_list->poly[i]->vlist->vtx[j]->x = posX;
+					poly_list->poly[i]->vlist->vtx[j]->y = posY;
+					poly_list->poly[i]->vlist->vtx[j]->z = posZ;
+					cellidx=vertex_self_avoidance(vesicle, poly_list->poly[i]->vlist->vtx[j]);
+					retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,poly_list->poly[i]->vlist->vtx[j]);
+					if(retval==TS_SUCCESS){
+						retval=cell_add_vertex(vesicle->clist->cell[cellidx],poly_list->poly[i]->vlist->vtx[j]);
+						break;
+					}
+					else{
+					//	printf("%d %d Cannot put the vertex here. Finding another position\n",i,j);
+						//randomly change the direction.
+						m=0;
+						//we must move first vertex into the vesicle if the normal is in or out of the vesicle if the normal is out
+						do{
+							costheta=2.0*drand48()-1.0;
+							sintheta=sqrt(1-pow(costheta,2));
+							phi=drand48()*2.0*M_PI;
+							if(j==0){
+						//for special cases, when we are on the edge of bipyramid  the distance od dmin_interspecies is not enough
+								posX=prevPosX+vesicle->dmax*sintheta*cos(phi);
+								posY=prevPosY+vesicle->dmax*sintheta*sin(phi);
+								posZ=prevPosZ+vesicle->dmax*costheta;
+
+							} else {
+								posX=prevPosX+vesicle->clist->dmin_interspecies*sintheta*cos(phi);
+								posY=prevPosY+vesicle->clist->dmin_interspecies*sintheta*sin(phi);
+								posZ=prevPosZ+vesicle->clist->dmin_interspecies*costheta;
+							}
+							m++;
+							if(m>1000) {
+								k=9999; //break also ot of the outer loop
+								printf("was here\n");
+								break;
+							}
+						}
+						while((xnorm*(poly_list->poly[i]->grafted_vtx->x-posX)+ynorm*(poly_list->poly[i]->grafted_vtx->y-posY)+znorm*(poly_list->poly[i]->grafted_vtx->z-posZ))>0.0 && j==0);
+					}
+					k++;
+					if(k>1000){
+						//undo changes to the cell
+						for(l=0;l<j;l++){
+							cellidx=vertex_self_avoidance(vesicle, poly_list->poly[i]->vlist->vtx[l]);
+							cell_remove_vertex(vesicle->clist->cell[cellidx],poly_list->poly[i]->vlist->vtx[l]);
+						}
+						return TS_FAIL;
+					}
+				}
+				prevPosX=posX;
+				prevPosY=posY;
+				prevPosZ=posZ;
+			}
+	printf("did it\n");
+	return TS_SUCCESS;
+
+}
+
+
+ts_poly_list *init_poly_list(ts_uint n_poly, ts_uint n_mono, ts_vertex_list *vlist, ts_vesicle *vesicle){
+	ts_poly_list *poly_list=(ts_poly_list *)calloc(1,sizeof(ts_poly_list));
+	poly_list->poly	= (ts_poly **)calloc(n_poly,sizeof(ts_poly *));
+	ts_uint i=0,j=0; //idx;
+	ts_uint gvtxi;
+	ts_bool retval;
+	ts_double dphi,dh;
+	cell_occupation(vesicle); //needed for evading the membrane
+	// Grafting polymers:
+	int tries=0;
+	if (vlist!=NULL){
+		if (n_poly > vlist->n){fatal("Number of polymers larger than numbero f vertices on a vesicle.",310);}
+		while(i<n_poly){
+			gvtxi = rand() % vlist->n;
+			if (vlist->vtx[gvtxi]->grafted_poly == NULL){
+				poly_list->poly[i] = init_poly(n_mono, vlist->vtx[gvtxi]);
+				retval=poly_initial_distribution(poly_list, i, vesicle);
+				if(retval==TS_FAIL){
+					ts_fprintf(stdout,"Found new potential grafting vertex %d for poly %d\n",gvtxi,i);
+					poly_free(poly_list->poly[i]);
+					tries++;
+				}
+				else {
+					tries=0;
+					i++;
+				}
+				if(tries>5000){
+					fatal("Cannot find space for inner polymeres",1001);
+				}
 			}
 		}
+	}
+	else
+	{
+		for(i=0;i<n_poly;i++){
+			poly_list->poly[i] = init_poly(n_mono, NULL);
+		}
+	}
+
+	poly_list->n = n_poly;
+
+	if (vlist!=NULL){
 	}
 	else
 	{
@@ -149,7 +235,7 @@
 		}
 	}
 */
-
+	
 	return poly_list;
 }
 
diff --git a/src/tape b/src/tape
index a7baa73..93e42d7 100644
--- a/src/tape
+++ b/src/tape
@@ -25,9 +25,9 @@
 
 ####### Polymer (brush) definitions ###########
 # npoly is a number of polymers attached to npoly distinct vertices on vesicle
-npoly=800
+npoly=400
 # nmono is a number of monomers in each polymer
-nmono=6
+nmono=10
 # Spring constant between monomers of the polymer
 k_spring=800
 #set to 1 if half of the polymeres are inside the vesicle
@@ -59,7 +59,7 @@
 mcsweeps=200
 #how many initial mcsweeps*inititer MC sweeps before recording to disk?
 #2
-inititer=10
+inititer=1
 #how many records do you want on the disk iteration are there in a run?
 #10000
 iterations=100

--
Gitblit v1.9.3