From 96b0a0ac783d35fbf4e764fe8e610ddc2a79af27 Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Wed, 04 Jan 2017 22:03:25 +0000 Subject: [PATCH] Added even more readability... --- src/bondflip.c | 151 +++++++++++++++++++++++++++++++++++++------------ 1 files changed, 113 insertions(+), 38 deletions(-) diff --git a/src/bondflip.c b/src/bondflip.c index 31836e8..66c40f7 100644 --- a/src/bondflip.c +++ b/src/bondflip.c @@ -1,3 +1,4 @@ +/* vim: set ts=4 sts=4 sw=4 noet : */ #include<stdlib.h> #include<math.h> #include "general.h" @@ -31,7 +32,7 @@ ts_vertex *k=bond->vtx2; ts_uint nei,neip,neim; ts_uint i,j; - ts_double oldenergy, delta_energy, dvol=0.0; + ts_double oldenergy, delta_energy, dvol=0.0, darea=0.0; ts_triangle *lm=NULL,*lp=NULL, *lp1=NULL, *lm2=NULL; ts_vertex *kp,*km; @@ -167,12 +168,17 @@ oldenergy+=kp->xk* kp->energy; oldenergy+=km->xk* km->energy; oldenergy+=it->xk* it->energy; + oldenergy+=bond->energy; /* attraction with neighboring vertices, that have spontaneous curvature */ //Neigbours of k, it, km, kp don't change its energy. - if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch==1){dvol = -lm->volume - lp->volume;} + if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch>0){dvol = -lm->volume - lp->volume;} + if(vesicle->tape->constareaswitch==2){darea=-lm->area-lp->area;} +/* vesicle_volume(vesicle); + fprintf(stderr,"Volume in the beginning=%1.16e\n", vesicle->volume); +*/ /* fix data structure for flipped bond */ - ts_flip_bond(k,it,km,kp, bond,lm, lp, lm2, lp1); + ts_flip_bond(k,it,km,kp, bond,lm, lp, lm2, lp1, vesicle->tape->w); /* Calculating the new energy */ @@ -181,46 +187,105 @@ delta_energy+=kp->xk* kp->energy; delta_energy+=km->xk* km->energy; delta_energy+=it->xk* it->energy; + delta_energy+=bond->energy; /* attraction with neighboring vertices, that have spontaneous curvature */ //Neigbours of k, it, km, kp don't change its energy. delta_energy-=oldenergy; - if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch==1){ + if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch>0){ dvol = dvol + lm->volume + lp->volume; if(vesicle->pswitch==1) delta_energy-= vesicle->pressure*dvol; } + if(vesicle->tape->constareaswitch==2){ + darea=darea+lm->area+lp->area; +/*check whether the dvol is gt than epsvol */ + if(fabs(vesicle->area+darea-A0)>epsarea){ + //restore old state. + /* restoration procedure copied from few lines below */ + for(i=0;i<4;i++){ + // fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); + free(orig_vtx[i]->neigh); + free(orig_vtx[i]->tristar); + free(orig_vtx[i]->bond); + free(orig_tria[i]->neigh); + memcpy((void *)orig_vtx[i],(void *)bck_vtx[i],sizeof(ts_vertex)); + memcpy((void *)orig_tria[i],(void *)bck_tria[i],sizeof(ts_triangle)); + // fprintf(stderr,"Restored vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); + /* level 2 pointers are redirected*/ + } + memcpy(bond,bck_bond,sizeof(ts_bond)); + for(i=0;i<4;i++){ + free(bck_vtx[i]); + free(bck_tria[i]); + /* fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d =",i, orig_vtx[i]->neigh_no ); + for(j=0;j<orig_vtx[i]->neigh_no;j++) fprintf(stderr," %d", orig_vtx[i]->neigh[j]->idx); + fprintf(stderr,"\n"); */ + } + free(bck_bond); + return TS_FAIL; - retval=TS_SUCCESS; + } + } + + + + + if(vesicle->tape->constvolswitch == 2){ + /*check whether the dvol is gt than epsvol */ + if(fabs(vesicle->volume+dvol-V0)>epsvol){ + //restore old state. + /* restoration procedure copied from few lines below */ + for(i=0;i<4;i++){ + // fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); + free(orig_vtx[i]->neigh); + free(orig_vtx[i]->tristar); + free(orig_vtx[i]->bond); + free(orig_tria[i]->neigh); + memcpy((void *)orig_vtx[i],(void *)bck_vtx[i],sizeof(ts_vertex)); + memcpy((void *)orig_tria[i],(void *)bck_tria[i],sizeof(ts_triangle)); + // fprintf(stderr,"Restored vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); + /* level 2 pointers are redirected*/ + } + memcpy(bond,bck_bond,sizeof(ts_bond)); + for(i=0;i<4;i++){ + free(bck_vtx[i]); + free(bck_tria[i]); + /* fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d =",i, orig_vtx[i]->neigh_no ); + for(j=0;j<orig_vtx[i]->neigh_no;j++) fprintf(stderr," %d", orig_vtx[i]->neigh[j]->idx); + fprintf(stderr,"\n"); */ + } + free(bck_bond); + return TS_FAIL; + + } + + } else if(vesicle->tape->constvolswitch == 1){ retval=constvolume(vesicle, it, -dvol, &delta_energy_cv, &constvol_vtx_moved,&constvol_vtx_backup); if(retval==TS_FAIL){ /* restoration procedure copied from few lines below */ - for(i=0;i<4;i++){ -// fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); - free(orig_vtx[i]->neigh); - free(orig_vtx[i]->tristar); - free(orig_vtx[i]->bond); - free(orig_tria[i]->neigh); - memcpy((void *)orig_vtx[i],(void *)bck_vtx[i],sizeof(ts_vertex)); - memcpy((void *)orig_tria[i],(void *)bck_tria[i],sizeof(ts_triangle)); -// fprintf(stderr,"Restored vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); - /* level 2 pointers are redirected*/ - } - memcpy(bond,bck_bond,sizeof(ts_bond)); - - for(i=0;i<4;i++){ - free(bck_vtx[i]); - free(bck_tria[i]); -/* fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d =",i, orig_vtx[i]->neigh_no ); - for(j=0;j<orig_vtx[i]->neigh_no;j++) fprintf(stderr," %d", orig_vtx[i]->neigh[j]->idx); - fprintf(stderr,"\n"); */ - } - - free(bck_bond); - return TS_FAIL; - + for(i=0;i<4;i++){ + // fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); + free(orig_vtx[i]->neigh); + free(orig_vtx[i]->tristar); + free(orig_vtx[i]->bond); + free(orig_tria[i]->neigh); + memcpy((void *)orig_vtx[i],(void *)bck_vtx[i],sizeof(ts_vertex)); + memcpy((void *)orig_tria[i],(void *)bck_tria[i],sizeof(ts_triangle)); + // fprintf(stderr,"Restored vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); + /* level 2 pointers are redirected*/ + } + memcpy(bond,bck_bond,sizeof(ts_bond)); + for(i=0;i<4;i++){ + free(bck_vtx[i]); + free(bck_tria[i]); + /* fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d =",i, orig_vtx[i]->neigh_no ); + for(j=0;j<orig_vtx[i]->neigh_no;j++) fprintf(stderr," %d", orig_vtx[i]->neigh[j]->idx); + fprintf(stderr,"\n"); */ + } + free(bck_bond); + return TS_FAIL; } - - delta_energy+=delta_energy_cv; + delta_energy+=delta_energy_cv; } @@ -239,6 +304,9 @@ //not accepted, reverting changes //restore all backups // fprintf(stderr,"Restoring!!!\n"); + if(vesicle->tape->constvolswitch == 1){ + constvolumerestore(constvol_vtx_moved,constvol_vtx_backup); + } for(i=0;i<4;i++){ // fprintf(stderr,"Restoring vtx neigh[%d] with neighbours %d\n",i, orig_vtx[i]->neigh_no ); @@ -263,10 +331,9 @@ free(bck_bond); - if(vesicle->tape->constvolswitch == 1){ - constvolumerestore(constvol_vtx_moved,constvol_vtx_backup); - } // fprintf(stderr,"Restoration complete!!!\n"); +// vesicle_volume(vesicle); +// fprintf(stderr,"Volume after fail=%1.16e\n", vesicle->volume); return TS_FAIL; } @@ -274,6 +341,14 @@ /* IF BONDFLIP ACCEPTED, THEN RETURN SUCCESS! */ // fprintf(stderr,"SUCCESS!!!\n"); + if(vesicle->tape->constvolswitch == 2){ + vesicle->volume+=dvol; + } else if(vesicle->tape->constvolswitch == 1){ + constvolumeaccept(vesicle,constvol_vtx_moved,constvol_vtx_backup); + } + if(vesicle->tape->constareaswitch==2){ + vesicle->area+=darea; + } // delete all backups for(i=0;i<4;i++){ free(bck_vtx[i]->neigh); @@ -288,16 +363,15 @@ */ } free(bck_bond); - if(vesicle->tape->constvolswitch == 1){ - constvolumeaccept(vesicle,constvol_vtx_moved,constvol_vtx_backup); - } +// vesicle_volume(vesicle); +// fprintf(stderr,"Volume after success=%1.16e\n", vesicle->volume); return TS_SUCCESS; } ts_bool ts_flip_bond(ts_vertex *k,ts_vertex *it,ts_vertex *km, ts_vertex *kp, -ts_bond *bond, ts_triangle *lm, ts_triangle *lp, ts_triangle *lm2, ts_triangle *lp1){ +ts_bond *bond, ts_triangle *lm, ts_triangle *lp, ts_triangle *lm2, ts_triangle *lp1, ts_double w_energy){ ts_uint i; //lmidx, lpidx; if(k==NULL || it==NULL || km==NULL || kp==NULL){ @@ -363,6 +437,7 @@ energy_vertex(kp); energy_vertex(km); energy_vertex(it); + attraction_bond_energy(bond, w_energy); // END modifications to data structure! return TS_SUCCESS; } -- Gitblit v1.9.3