From a63f1719d2c7fd2c69accc0eb3eb038af50e555e Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Fri, 13 Jul 2012 13:10:01 +0000 Subject: [PATCH] Fixed frame to make changes in cell. --- src/vertexmove.c | 111 ++++++++++++++++++++++++++++++------------------------- 1 files changed, 61 insertions(+), 50 deletions(-) diff --git a/src/vertexmove.c b/src/vertexmove.c index 243463f..f38e3c7 100644 --- a/src/vertexmove.c +++ b/src/vertexmove.c @@ -11,60 +11,72 @@ //#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){ ts_uint i; ts_double dist; - ts_vertex *tvtx=(ts_vertex *)malloc(sizeof(ts_vertex)); - tvtx->data=init_vertex_data(); ts_bool retval; ts_uint cellidx; - ts_double xold,yold,zold; ts_double delta_energy,oenergy; - ts_vertex *ovtx; + //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->data->x=vtx->data->x+vesicle->stepsize*(2.0*rn[0]-1.0); - tvtx->data->y=vtx->data->y+vesicle->stepsize*(2.0*rn[1]-1.0); - tvtx->data->z=vtx->data->z+vesicle->stepsize*(2.0*rn[2]-1.0); - //check we if some length to neighbours are too much - for(i=0;i<vtx->data->neigh_no;i++){ - dist=vtx_distance_sq(tvtx,vtx->data->neigh[i]); - if(dist<1.0 || dist>vesicle->dmax) return TS_FAIL; + //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); + //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; + } } -fprintf(stderr,"Was here!\n"); //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=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->data->x; - yold=vtx->data->y; - zold=vtx->data->z; - ovtx=malloc(sizeof(ts_vertex)); - vtx_copy(ovtx,vtx); - vtx->data->x=tvtx->data->x; - vtx->data->y=tvtx->data->y; - vtx->data->z=tvtx->data->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->data->tristar_no;i++) triangle_normal_vector(vtx->data->tristar[i]); - //energy and curvature + for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); + oenergy=vtx->energy; energy_vertex(vtx); - delta_energy=vtx->data->xk*(vtx->data->energy - ovtx->data->energy); + delta_energy=vtx->xk*(vtx->energy - oenergy); //the same is done for neighbouring vertices - for(i=0;i<vtx->data->neigh_no;i++){ - oenergy=vtx->data->neigh[i]->data->energy; - energy_vertex(vtx->data->neigh[i]); - delta_energy+=vtx->data->neigh[i]->data->xk*(vtx->data->neigh[i]->data->energy-oenergy); + 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); } - fprintf(stderr, "DE=%f\n",delta_energy); +// fprintf(stderr, "DE=%f\n",delta_energy); //MONTE CARLOOOOOOOO if(delta_energy>=0){ #ifdef TS_DOUBLE_DOUBLE @@ -78,28 +90,27 @@ #endif { //not accepted, reverting changes - vtx->data->x=xold; - vtx->data->y=yold; - vtx->data->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->data->tristar_no;i++) triangle_normal_vector(vtx->data->tristar[i]); - //energy and curvature - energy_vertex(vtx); - //the same is done for neighbouring vertices - for(i=0;i<vtx->data->neigh_no;i++) energy_vertex(vtx->data->neigh[i]); - free(ovtx->data->bond_length); - free(ovtx->data->bond_length_dual); - free(ovtx); + 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! - - free(ovtx->data->bond_length); - free(ovtx->data->bond_length_dual); - free(ovtx); return TS_SUCCESS; } -- Gitblit v1.9.3