From ded11ede7e9eb11c216d35d61ab9331e697f47bc Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Sat, 07 Dec 2013 09:09:37 +0000 Subject: [PATCH] commiting before checking out onother branch --- src/vertexmove.c | 102 +++++++++++++++++++++++++++++---------------------- 1 files changed, 58 insertions(+), 44 deletions(-) diff --git a/src/vertexmove.c b/src/vertexmove.c index b575346..7047450 100644 --- a/src/vertexmove.c +++ b/src/vertexmove.c @@ -20,55 +20,73 @@ ts_bool retval; ts_uint cellidx; ts_double delta_energy,oenergy; + ts_double costheta,sintheta,phi,r; //This will hold all the information of vtx and its neighbours - ts_vertex **backupvtx=(ts_vertex **)calloc(vtx->neigh_no+1,sizeof(ts_vertex *)); - backupvtx[0]=(ts_vertex *)malloc(sizeof(ts_vertex)); - backupvtx[0]=(ts_vertex *)memcpy((void *)backupvtx[0],(void *)vtx,sizeof(ts_vertex)); + ts_vertex backupvtx[20]; + memcpy((void *)&backupvtx[0],(void *)vtx,sizeof(ts_vertex)); + + //Some stupid tests for debugging cell occupation! +/* cellidx=vertex_self_avoidance(vesicle, vtx); + if(vesicle->clist->cell[cellidx]==vtx->cell){ + fprintf(stderr,"Idx match!\n"); + } else { + fprintf(stderr,"***** Idx don't match!\n"); + fatal("ENding.",1); + } +*/ + //temporarly moving the vertex - vtx->x=vtx->x+vesicle->stepsize*(2.0*rn[0]-1.0); - vtx->y=vtx->y+vesicle->stepsize*(2.0*rn[1]-1.0); - vtx->z=vtx->z+vesicle->stepsize*(2.0*rn[2]-1.0); - //check we if some length to neighbours are too much - for(i=0;i<vtx->neigh_no;i++){ - dist=vtx_distance_sq(vtx,vtx->neigh[i]); +// vtx->x=vtx->x+vesicle->stepsize*(2.0*rn[0]-1.0); +// 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: + r=vesicle->stepsize*rn[0]; + phi=rn[1]*2*M_PI; + costheta=2*rn[2]-1; + sintheta=sqrt(1-pow(costheta,2)); + vtx->x=vtx->x+r*sintheta*cos(phi); + vtx->y=vtx->y+r*sintheta*sin(phi); + vtx->z=vtx->z+r*costheta; + + + //distance with neighbours check + for(i=0;i<vtx->neigh->n;i++){ + dist=vtx_distance_sq(vtx,vtx->neigh->vtx[i]); if(dist<1.0 || dist>vesicle->dmax) { - vtx=memcpy((void *)vtx,(void *)backupvtx[0],sizeof(ts_vertex)); - free(backupvtx[0]); - free(backupvtx); -// fprintf(stderr,"Fail 1, dist=%f, vesicle->dmax=%f\n", dist, vesicle->dmax); + 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); //check occupation number - retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,backupvtx[0],vtx); + retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,vtx); + if(retval==TS_FAIL){ - vtx=memcpy((void *)vtx,(void *)backupvtx[0],sizeof(ts_vertex)); - free(backupvtx[0]); - free(backupvtx); -// fprintf(stderr,"Fail 2\n"); + vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); return TS_FAIL; } //if all the tests are successful, then energy for vtx and neighbours is calculated - for(i=0;i<vtx->neigh_no;i++){ - backupvtx[i+1]=(ts_vertex *)malloc(sizeof(ts_vertex)); - backupvtx[i+1]=memcpy((void *)backupvtx[i+1],(void *)vtx->neigh[i],sizeof(ts_vertex)); + for(i=0;i<vtx->neigh->n;i++){ + memcpy((void *)&backupvtx[i+1],(void *)vtx->neigh->vtx[i],sizeof(ts_vertex)); } + + delta_energy=0; //update the normals of triangles that share bead i. for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); - //energy and curvature + oenergy=vtx->energy; energy_vertex(vtx); - delta_energy=vtx->xk*(vtx->energy - backupvtx[0]->energy); + delta_energy=vtx->xk*(vtx->energy - oenergy); //the same is done for neighbouring vertices - for(i=0;i<vtx->neigh_no;i++){ - oenergy=vtx->neigh[i]->energy; - energy_vertex(vtx->neigh[i]); - delta_energy+=vtx->neigh[i]->xk*(vtx->neigh[i]->energy-oenergy); + for(i=0;i<vtx->neigh->n;i++){ + oenergy=vtx->neigh->vtx[i]->energy; + energy_vertex(vtx->neigh->vtx[i]); + delta_energy+=vtx->neigh->vtx[i]->xk*(vtx->neigh->vtx[i]->energy-oenergy); } // fprintf(stderr, "DE=%f\n",delta_energy); //MONTE CARLOOOOOOOO @@ -84,31 +102,27 @@ #endif { //not accepted, reverting changes - vtx=memcpy((void *)vtx,(void *)backupvtx[0],sizeof(ts_vertex)); - free(backupvtx[0]); - for(i=0;i<vtx->neigh_no;i++){ - vtx->neigh[i]=memcpy((void *)vtx->neigh[i],(void *)backupvtx[i+1],sizeof(ts_vertex)); - free(backupvtx[i+1]); + vtx=memcpy((void *)vtx,(void *)&backupvtx[0],sizeof(ts_vertex)); + for(i=0;i<vtx->neigh->n;i++){ + vtx->neigh->vtx[i]=memcpy((void *)vtx->neigh->vtx[i],(void *)&backupvtx[i+1],sizeof(ts_vertex)); } - free(backupvtx); -// fprintf(stderr,"Reverted\n"); //update the normals of triangles that share bead i. - for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); + for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); return TS_FAIL; } } - //END MONTE CARLOOOOOOO - - //TODO: change cell occupation if necessary! -// fprintf(stderr,"Success!!\n"); - free(backupvtx[0]); - for(i=0;i<vtx->neigh_no;i++){ - free(backupvtx[i+1]); + +// oldcellidx=vertex_self_avoidance(vesicle, &backupvtx[0]); + 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); + } - free(backupvtx); -// fprintf(stderr,"Accepted\n"); +// if(oldcellidx); + //END MONTE CARLOOOOOOO return TS_SUCCESS; } -- Gitblit v1.9.3