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