From 311301bef2f4762781817a7553b56d827d2d8d7e Mon Sep 17 00:00:00 2001 From: Miha <mihaf@hapi.(none)> Date: Fri, 29 Nov 2013 08:44:43 +0000 Subject: [PATCH] First test of new branch. Nothing done yet. --- src/vertexmove.c | 106 ++++++++++++++++++++++++++++------------------------ 1 files changed, 57 insertions(+), 49 deletions(-) diff --git a/src/vertexmove.c b/src/vertexmove.c index 7edbf04..a79b2f0 100644 --- a/src/vertexmove.c +++ b/src/vertexmove.c @@ -11,6 +11,7 @@ //#include "io.h" #include<stdio.h> #include "vertexmove.h" +#include <string.h> ts_bool single_verticle_timestep(ts_vesicle *vesicle,ts_vertex *vtx,ts_double *rn){ @@ -18,59 +19,69 @@ ts_double dist; ts_bool retval; ts_uint cellidx; - ts_double xold,yold,zold; ts_double delta_energy,oenergy; - ts_vertex *ovtx; - ts_vertex *tvtx=(ts_vertex *)calloc(1,sizeof(ts_vertex)); ts_double costheta,sintheta,phi,r; + //This will hold all the information of vtx and its neighbours + ts_vertex backupvtx[20]; + memcpy((void *)&backupvtx[0],(void *)vtx,sizeof(ts_vertex)); - //randomly we move the temporary vertex -// tvtx->x=vtx->x+vesicle->stepsize*(2.0*rn[0]-1.0); -// tvtx->y=vtx->y+vesicle->stepsize*(2.0*rn[1]-1.0); -// tvtx->z=vtx->z+vesicle->stepsize*(2.0*rn[2]-1.0); + //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); + //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)); - tvtx->x=vtx->x+r*sintheta*cos(phi); - tvtx->y=vtx->y+r*sintheta*sin(phi); - tvtx->z=vtx->z+r*costheta; - //check we if some length to neighbours are too much + 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_no;i++){ - dist=vtx_distance_sq(tvtx,vtx->neigh[i]); + dist=vtx_distance_sq(vtx,vtx->neigh[i]); if(dist<1.0 || dist>vesicle->dmax) { - vtx_free(tvtx); -// 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, tvtx); + cellidx=vertex_self_avoidance(vesicle, vtx); //check occupation number - retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,vtx,tvtx); + retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,vtx); + if(retval==TS_FAIL){ - vtx_free(tvtx); -// 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 we update the vertex position - xold=vtx->x; - yold=vtx->y; - zold=vtx->z; - ovtx=malloc(sizeof(ts_vertex)); - vtx_copy(ovtx,vtx); - vtx->x=tvtx->x; - vtx->y=tvtx->y; - vtx->z=tvtx->z; + + + //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)); + } + + 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 - ovtx->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; @@ -91,30 +102,27 @@ #endif { //not accepted, reverting changes - vtx->x=xold; - vtx->y=yold; - vtx->z=zold; + 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]); - //energy and curvature - energy_vertex(vtx); - //the same is done for neighbouring vertices - for(i=0;i<vtx->neigh_no;i++) energy_vertex(vtx->neigh[i]); - free(ovtx->bond_length); - free(ovtx->bond_length_dual); - free(ovtx); - vtx_free(tvtx); + for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); + return TS_FAIL; } } + +// 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); + + } +// if(oldcellidx); //END MONTE CARLOOOOOOO - - //TODO: change cell occupation if necessary! -// fprintf(stderr,"Success!!\n"); - free(ovtx->bond_length); - free(ovtx->bond_length_dual); - free(ovtx); - vtx_free(tvtx); return TS_SUCCESS; } -- Gitblit v1.9.3