From aafd67f2d592799167864698ea07eb8aadcc86a7 Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Wed, 25 May 2016 18:37:21 +0000 Subject: [PATCH] Merge branch 'master' of bitbucket.org:samop/trisurf-ng with work done on removing the memory leaks. --- src/vertexmove.c | 74 ++++++++++++++++++++++++++++++++----- 1 files changed, 64 insertions(+), 10 deletions(-) diff --git a/src/vertexmove.c b/src/vertexmove.c index cecb73e..4caa86f 100644 --- a/src/vertexmove.c +++ b/src/vertexmove.c @@ -1,3 +1,4 @@ +/* vim: set ts=4 sts=4 sw=4 noet : */ #include<stdlib.h> #include<math.h> #include "general.h" @@ -20,7 +21,7 @@ ts_double dist; ts_bool retval; ts_uint cellidx; - ts_double delta_energy, delta_energy_cv,oenergy,dvol=0.0; + ts_double delta_energy, delta_energy_cv,oenergy,dvol=0.0, darea=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; @@ -41,7 +42,7 @@ // vtx->y=vtx->y+vesicle->stepsize*(2.0*rn[1]-1.0); // vtx->z=vtx->z+vesicle->stepsize*(2.0*rn[2]-1.0); - //random move in a sphere with radius stepsize: +//random move in a sphere with radius stepsize: r=vesicle->stepsize*rn[0]; phi=rn[1]*2*M_PI; costheta=2*rn[2]-1; @@ -51,12 +52,12 @@ vtx->z=vtx->z+r*costheta; - //distance with neighbours check +//distance with neighbours check for(i=0;i<vtx->neigh_no;i++){ dist=vtx_distance_sq(vtx,vtx->neigh[i]); if(dist<1.0 || dist>vesicle->dmax) { - vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); - return TS_FAIL; + vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); + return TS_FAIL; } } @@ -71,10 +72,19 @@ // TODO: Maybe faster if checks only nucleus-neighboring cells // Nucleus penetration check: +if(vesicle->R_nucleus>0.0){ if (vtx->x*vtx->x + vtx->y*vtx->y + vtx->z*vtx->z < vesicle->R_nucleus){ vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); return TS_FAIL; } +} else if(vesicle->R_nucleusX>0.0){ +// fprintf(stderr,"DEBUG, (Rx, Ry,Rz)^2=(%f,%f,%f)\n",vesicle->R_nucleusX, vesicle->R_nucleusY, vesicle->R_nucleusZ); + if ((vtx->x*vtx->x)/vesicle->R_nucleusX + vtx->y*vtx->y/vesicle->R_nucleusY + (vtx->z*vtx->z)/vesicle->R_nucleusZ < 1.0){ + vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); + return TS_FAIL; + } + +} //self avoidance check with distant vertices cellidx=vertex_self_avoidance(vesicle, vtx); @@ -87,14 +97,19 @@ } - //if all the tests are successful, then energy for vtx and neighbours is calculated +//if all the tests are successful, then energy for vtx and neighbours is calculated for(i=0;i<vtx->neigh_no;i++){ memcpy((void *)&backupvtx[i+1],(void *)vtx->neigh[i],sizeof(ts_vertex)); } - if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch==1){ + if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch>0){ for(i=0;i<vtx->tristar_no;i++) dvol-=vtx->tristar[i]->volume; - }; + } + + if(vesicle->tape->constareaswitch==2){ + for(i=0;i<vtx->tristar_no;i++) darea-=vtx->tristar[i]->area; + + } delta_energy=0; @@ -113,11 +128,43 @@ delta_energy+=vtx->neigh[i]->xk*(vtx->neigh[i]->energy-oenergy); } - if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch == 1){ + if(vesicle->pswitch == 1 || vesicle->tape->constvolswitch >0){ for(i=0;i<vtx->tristar_no;i++) dvol+=vtx->tristar[i]->volume; if(vesicle->pswitch==1) delta_energy-=vesicle->pressure*dvol; }; + 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; + } + + + } + + if(vesicle->tape->constvolswitch==2){ + /*check whether the dvol is gt than epsvol */ + //fprintf(stderr,"DVOL=%1.16e\n",dvol); + if(fabs(vesicle->volume+dvol-V0)>epsvol){ + //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; + } + + } else // vesicle_volume(vesicle); // fprintf(stderr,"Volume before=%1.16e\n", vesicle->volume); if(vesicle->tape->constvolswitch == 1){ @@ -150,7 +197,7 @@ //MONTE CARLOOOOOOOO if(delta_energy>=0){ #ifdef TS_DOUBLE_DOUBLE - if(exp(-delta_energy)< drand48() ) + if(exp(-delta_energy)< drand48()) #endif #ifdef TS_DOUBLE_FLOAT if(expf(-delta_energy)< (ts_float)drand48()) @@ -189,9 +236,16 @@ } + 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; + } // if(oldcellidx); //END MONTE CARLOOOOOOO // vesicle_volume(vesicle); -- Gitblit v1.9.3