From 3b39021cfdf157d35ce3475fb33b2c7a35bce676 Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Wed, 27 Jun 2012 14:09:49 +0000 Subject: [PATCH] Debugged the spherical harmonics --- src/vertexmove.c | 70 +++++++++++++++++++---------------- 1 files changed, 38 insertions(+), 32 deletions(-) diff --git a/src/vertexmove.c b/src/vertexmove.c index 243463f..c62aeda 100644 --- a/src/vertexmove.c +++ b/src/vertexmove.c @@ -16,55 +16,59 @@ *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; + ts_vertex *tvtx=(ts_vertex *)calloc(1,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); + 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); //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; + for(i=0;i<vtx->neigh_no;i++){ + dist=vtx_distance_sq(tvtx,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); + return TS_FAIL; + } } -fprintf(stderr,"Was here!\n"); //self avoidance check with distant vertices cellidx=vertex_self_avoidance(vesicle, tvtx); //check occupation number retval=cell_occupation_number_and_internal_proximity(vesicle->clist,cellidx,vtx,tvtx); if(retval==TS_FAIL){ + vtx_free(tvtx); +// fprintf(stderr,"Fail 2\n"); 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; + xold=vtx->x; + yold=vtx->y; + zold=vtx->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; + vtx->x=tvtx->x; + vtx->y=tvtx->y; + vtx->z=tvtx->z; 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]); + for(i=0;i<vtx->tristar_no;i++) triangle_normal_vector(vtx->tristar[i]); //energy and curvature energy_vertex(vtx); - delta_energy=vtx->data->xk*(vtx->data->energy - ovtx->data->energy); + delta_energy=vtx->xk*(vtx->energy - ovtx->energy); //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 +82,30 @@ #endif { //not accepted, reverting changes - vtx->data->x=xold; - vtx->data->y=yold; - vtx->data->z=zold; + vtx->x=xold; + vtx->y=yold; + vtx->z=zold; //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]); + 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->data->neigh_no;i++) energy_vertex(vtx->data->neigh[i]); - free(ovtx->data->bond_length); - free(ovtx->data->bond_length_dual); + 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); return TS_FAIL; } } //END MONTE CARLOOOOOOO //TODO: change cell occupation if necessary! - - free(ovtx->data->bond_length); - free(ovtx->data->bond_length_dual); +// 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