From 06f70cc6f065746cf79c9792228356e3f015f4ec Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Fri, 08 Mar 2019 20:13:00 +0000 Subject: [PATCH] Constant area plugin done -- fix --- src/vertexmove.c | 83 +++++++++++------------------------------ 1 files changed, 23 insertions(+), 60 deletions(-) diff --git a/src/vertexmove.c b/src/vertexmove.c index 37508e9..5a32955 100644 --- a/src/vertexmove.c +++ b/src/vertexmove.c @@ -20,7 +20,7 @@ ts_uint i; ts_bool retval; ts_uint cellidx; - ts_double delta_energy, oenergy,dvol=0.0, darea=0.0, dstretchenergy=0.0; + ts_double delta_energy, oenergy, darea=0.0, dstretchenergy=0.0; ts_double costheta,sintheta,phi,r; //This will hold all the information of vtx and its neighbours ts_vertex backupvtx[20]; // *constvol_vtx_moved=NULL, *constvol_vtx_backup=NULL; @@ -60,22 +60,15 @@ vesicle->plist->pointer->plugin->function->vm_energy_before_prepare(vesicle, vtx); vesicle->plist->pointer=vesicle->plist->pointer->next; } +/* End of vm_energy_before_prepare() */ - if(vesicle->tape->constareaswitch==2){ - for(i=0;i<vtx->tristar_no;i++) darea-=vtx->tristar[i]->area; - - } //stretching energy 1 of 3 if(vesicle->tape->stretchswitch==1){ for(i=0;i<vtx->tristar_no;i++) dstretchenergy-=vtx->tristar[i]->energy; } delta_energy=0; - -// vesicle_volume(vesicle); -// fprintf(stderr,"Volume in the beginning=%1.16e\n", vesicle->volume); - - //update the normals of triangles that share bead i. +//update the normals of triangles that share bead i. for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); oenergy=vtx->energy; energy_vertex(vtx); @@ -96,24 +89,6 @@ vesicle->plist->pointer=vesicle->plist->pointer->next; } - - - if(vesicle->tape->constareaswitch==2){ - /* check whether the darea is gt epsarea */ - for(i=0;i<vtx->tristar_no;i++) darea+=vtx->tristar[i]->area; - if(fabs(vesicle->area+darea-A0)>epsarea){ - //restore old state. - vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); - for(i=0;i<vtx->neigh_no;i++){ - vtx->neigh[i]=memcpy((void *)vtx->neigh[i],(void *)&backupvtx[i+1],sizeof(ts_vertex)); - } - for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); - //fprintf(stderr,"fajlam!\n"); - return TS_FAIL; - } - - - } /* Vertices with spontaneous curvature may have spontaneous force perpendicular to the surface of the vesicle. additional delta energy is calculated in this function */ delta_energy+=direct_force_energy(vesicle,vtx,backupvtx); @@ -185,15 +160,15 @@ if(expl(-delta_energy)< (ts_ldouble)drand48()) #endif { - //not accepted, reverting changes - // fprintf(stderr,"MC failed\n"); +/*************************************************** MC step rejected **************************************************************/ vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); for(i=0;i<vtx->neigh_no;i++){ vtx->neigh[i]=memcpy((void *)vtx->neigh[i],(void *)&backupvtx[i+1],sizeof(ts_vertex)); } - //update the normals of triangles that share bead i. - for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); + //update the normals of triangles that share bead i. + for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); + //stretching energy 3 of 3 if(vesicle->tape->stretchswitch==1){ for(i=0;i<vtx->tristar_no;i++){ @@ -201,37 +176,26 @@ } } -// fprintf(stderr, "before vtx(x,y,z)=%e,%e,%e\n",constvol_vtx_moved->x, constvol_vtx_moved->y, constvol_vtx_moved->z); -// if(vesicle->tape->constvolswitch == 1){ -// constvolumerestore(constvol_vtx_moved,constvol_vtx_backup); -// } -// fprintf(stderr, "after vtx(x,y,z)=%e,%e,%e\n",constvol_vtx_moved->x, constvol_vtx_moved->y, constvol_vtx_moved->z); -// vesicle_volume(vesicle); -// fprintf(stderr,"Volume after fail=%1.16e\n", vesicle->volume); + + +/* Entry point for plugin vm_before_montecarlo_constraint() function */ + vesicle->plist->pointer=vesicle->plist->chain->vm_new_state_rejected; + while(vesicle->plist->pointer!=NULL){ + vesicle->plist->pointer->plugin->function->vm_new_state_rejected(vesicle,vtx, &backupvtx[0]); + vesicle->plist->pointer=vesicle->plist->pointer->next; + } +/* End of vm_before_montecarlo_constraint() */ + + return TS_FAIL; } } - //accepted - // fprintf(stderr,"MC accepted\n"); -// oldcellidx=vertex_self_avoidance(vesicle, &backupvtx[0]); +/*************************************************** MC step accepted **************************************************************/ cellidx=vertex_self_avoidance(vesicle, vtx); if(vtx->cell!=vesicle->clist->cell[cellidx]){ retval=cell_add_vertex(vesicle->clist->cell[cellidx],vtx); -// if(retval==TS_SUCCESS) cell_remove_vertex(vesicle->clist->cell[oldcellidx],vtx); - if(retval==TS_SUCCESS) cell_remove_vertex(backupvtx[0].cell,vtx); - + if(retval==TS_SUCCESS) cell_remove_vertex(backupvtx[0].cell,vtx); } - - if(vesicle->tape->constvolswitch == 2){ - vesicle->volume+=dvol; - } - - - - if(vesicle->tape->constareaswitch==2){ - vesicle->area+=darea; - } - /* Entry point for plugin vm_before_montecarlo_constraint() function */ vesicle->plist->pointer=vesicle->plist->chain->vm_new_state_accepted; @@ -241,14 +205,13 @@ } /* End of vm_before_montecarlo_constraint() */ - - - - return TS_SUCCESS; } + + + ts_bool single_poly_vertex_move(ts_vesicle *vesicle,ts_poly *poly,ts_vertex *vtx,ts_double *rn){ ts_uint i; ts_bool retval; -- Gitblit v1.9.3