Trisurf Monte Carlo simulator
Samo Penic
2013-12-08 d33abe246221255be91bacda79f5852b8c4ed8fd
commit | author | age
d7639a 1 #include<stdlib.h>
SP 2 #include<math.h>
3 #include<string.h>
4 #include "general.h"
5 #include "vertex.h"
a10dd5 6 #include "bond.h"
d7639a 7 #include<stdio.h>
SP 8
73f967 9 ts_vertex_list *init_vertex_list(ts_uint N){    
d7639a 10     ts_int i;
bba37c 11     ts_vertex_list *vlist=(ts_vertex_list *)calloc(1,sizeof(ts_vertex_list));
SP 12  
d7639a 13     if(N==0){
bba37c 14             vlist->n=0;
SP 15         vlist->vtx=(ts_vertex **)calloc(1,sizeof(ts_vertex *)); /*Allocate one memory space for vlist anyway */
73f967 16         return vlist;
d7639a 17     }
73f967 18     
8f6a69 19     vlist->vtx=(ts_vertex **)calloc(N,sizeof(ts_vertex *));
SP 20     if(vlist->vtx==NULL)
73f967 21         fatal("Fatal error reserving memory space for vertex list! Could number of requsted vertices be too large?", 100);
a10dd5 22     for(i=0;i<N;i++) {
8f6a69 23         vlist->vtx[i]=(ts_vertex *)calloc(1,sizeof(ts_vertex));
a10dd5 24         vlist->vtx[i]->idx=i;
8ed8fd 25     vlist->vtx[i]->neigh=init_vertex_list(0);
bba37c 26         /* initialize Ylm for spherical hamonics DONE in sh.c */
SP 27         }
73f967 28     vlist->n=N;
ded11e 29     vlist->list_size=TS_VLIST_CHUNK; //TODO: can be buggy in some cases, when N>0 and we want to delete some vertices.
SP 30     return vlist;
d7639a 31 }
SP 32
64c113 33
SP 34 ts_bool vertex_list_add_vtx(ts_vertex_list *vlist, ts_vertex *vtx){
35
36 #ifdef DEBUG
8ed8fd 37     if(vtx==NULL || vlist==NULL)
SP 38         return TS_FAIL;
64c113 39 #endif
SP 40     if(vlist->list_size < vlist->n+1){
41         vlist->vtx=(ts_vertex **)realloc(vlist->vtx, (vlist->list_size+TS_VLIST_CHUNK)*sizeof(ts_vertex*));
42         if(vlist->vtx==NULL){
43             fatal("Error in vertex_list_add. Could not reallocate memory space.",9999);
44         }
45         vlist->list_size+=TS_VLIST_CHUNK;
46     }
bba37c 47     vlist->vtx[vlist->n]=vtx;
64c113 48     vlist->n++;
SP 49     return TS_SUCCESS;
50 }
51
52
53 /* Idea is to delete vertex by removing it from list. The emply place is then filled in by the last 
54 vertex in the list. List can then be resized by one -- unless resize is required when max list size - 
55 real list size > predefined value */
56 ts_bool vertex_list_remove_vtx(ts_vertex_list *vlist, ts_vertex *vtx){
57 #ifdef DEBUG
58     if(vtx==NULL || vlist==NULL) return TS_FAIL;
59 #endif
bba37c 60     ts_uint i;
64c113 61     for(i=0; i<vlist->n;i++){
SP 62         if(vlist->vtx[i]==vtx){
bba37c 63             vlist->vtx[i]=vlist->vtx[vlist->n-1];
64c113 64             vlist->n--;
SP 65             if(vlist->list_size-vlist->n > TS_VLIST_CHUNK){
66                 vlist->vtx=(ts_vertex **)realloc(vlist->vtx, (vlist->list_size-TS_VLIST_CHUNK)*sizeof(ts_vertex*));
67                 if(vlist->vtx==NULL){
68                     fatal("Error in vertex_list_add. Could not reallocate memory space.",9999);
69                 }
70                 vlist->list_size-=TS_VLIST_CHUNK;
71             }
72             return TS_SUCCESS;
73         }
74     }
75     return TS_FAIL;
9802f1 76 }
7958e9 77
9802f1 78
SP 79
a10dd5 80 ts_bool vtx_add_bond(ts_bond_list *blist,ts_vertex *vtx1,ts_vertex *vtx2){
SP 81     ts_bond *bond;
82     bond=bond_add(blist,vtx1,vtx2);
83     if(bond==NULL) return TS_FAIL;
8f6a69 84     vtx1->bond_no++;
30ee9c 85     vtx2->bond_no++;
3131dc 86    // vtx2->data->bond_no++;
a10dd5 87
8f6a69 88     vtx1->bond=(ts_bond **)realloc(vtx1->bond, vtx1->bond_no*sizeof(ts_bond *)); 
30ee9c 89     vtx2->bond=(ts_bond **)realloc(vtx2->bond, vtx2->bond_no*sizeof(ts_bond *)); 
3131dc 90    // vtx2->data->bond=(ts_bond **)realloc(vtx2->data->bond, vtx2->data->bond_no*sizeof(ts_bond *)); 
8f6a69 91     vtx1->bond[vtx1->bond_no-1]=bond;
30ee9c 92     vtx2->bond[vtx2->bond_no-1]=bond;
8f6a69 93    // vtx2->ata->bond[vtx2->data->bond_no-1]=bond;
a10dd5 94     return TS_SUCCESS;
SP 95 }
96
bc4e1e 97 // Add neighbour just in one direction (add vtx2 as a neigh of vtx1 and not another way around!
2f8ee7 98 ts_bool vtx_add_neighbour(ts_vertex *vtx1, ts_vertex *vtx2){
8ed8fd 99     ts_uint i;
SP 100     if(vtx1==NULL || vtx2==NULL) return TS_FAIL;
101     for(i=0;i<vtx1->neigh->n;i++){
102         if (vtx1->neigh->vtx[i]==vtx2){
103             return TS_FAIL;
104         }
105     }
106         vertex_list_add_vtx(vtx1->neigh, vtx2);
2f8ee7 107     return TS_SUCCESS;
SP 108 }
109
ded11e 110
2f8ee7 111 ts_bool vtx_add_cneighbour(ts_bond_list *blist, ts_vertex *vtx1, ts_vertex *vtx2){
a10dd5 112     ts_bool retval;
8ed8fd 113     retval=vtx_add_neighbour(vtx1, vtx2);
SP 114 //    retval=vertex_list_add_vtx(vtx1->neigh, vtx2);
115 //    retval=vertex_list_add_vtx(vtx2->neigh, vtx1);
116     if(retval==TS_SUCCESS){
117     //    fprintf(stderr,"Bond added\n");
118        retval=vtx_add_bond(blist,vtx1,vtx2); 
119         }
a10dd5 120     return retval;
SP 121 }
122
9802f1 123 /*TODO: write and optimize this urgently before use! */
SP 124 ts_bool vtx_remove_cneighbour(ts_bond_list *blist, ts_vertex *vtx1, ts_vertex
125 *vtx2){
126 //    ts_bool retval;
127 /* remove the bond */
128 //retval=vtx_remove_bond(blist,vtx1,vtx2);
129 /* remove the vertices */
130     return TS_SUCCESS;
131 }
a10dd5 132
d7639a 133
737714 134 ts_bool vtx_free(ts_vertex  *vtx){
8ed8fd 135     if(vtx->neigh!=NULL)  {
SP 136         if (vtx->neigh->vtx!=NULL) free(vtx->neigh->vtx);
137         free(vtx->neigh);
138     }
8f6a69 139     if(vtx->tristar!=NULL) free(vtx->tristar);
SP 140     if(vtx->bond!=NULL)    free(vtx->bond);
737714 141     free(vtx);
d7639a 142     return TS_SUCCESS;
SP 143 }
144
73f967 145 ts_bool vtx_list_free(ts_vertex_list *vlist){
SP 146     int i;
147     for(i=0;i<vlist->n;i++){
8f6a69 148         if(vlist->vtx[i]!=NULL) vtx_free(vlist->vtx[i]);
73f967 149     }
8f6a69 150     //free(*(vlist->vtx));
737714 151     free(vlist->vtx);
73f967 152     free(vlist);
SP 153     return TS_SUCCESS;
154 }
d7639a 155
bb77ca 156 inline ts_double vtx_distance_sq(ts_vertex *vtx1, ts_vertex *vtx2){
d7639a 157     ts_double dist;
SP 158 #ifdef TS_DOUBLE_DOUBLE
8f6a69 159     dist=pow(vtx1->x-vtx2->x,2) + pow(vtx1->y-vtx2->y,2) + pow(vtx1->z-vtx2->z,2);
d7639a 160 #endif
SP 161 #ifdef TS_DOUBLE_LONGDOUBLE
8f6a69 162     dist=powl(vtx1->x-vtx2->x,2) + powl(vtx1->y-vtx2->y,2) + powl(vtx1->z-vtx2->z,2);
d7639a 163 #endif
SP 164 #ifdef TS_DOUBLE_FLOAT
8f6a69 165     dist=powf(vtx1->x-vtx2->x,2) + powf(vtx1->y-vtx2->y,2) + powf(vtx1->z-vtx2->z,2);
d7639a 166 #endif
SP 167     return(dist);
168 }
bb77ca 169
7958e9 170
SP 171
172 ts_bool vtx_set_global_values(ts_vesicle *vesicle){ 
173     ts_double xk=vesicle->bending_rigidity;
174     ts_uint i; 
175     for(i=0;i<vesicle->vlist->n;i++){
8f6a69 176         vesicle->vlist->vtx[i]->xk=xk;
7958e9 177     }
SP 178     return TS_SUCCESS;
179 }
180
181 inline ts_double vtx_direct(ts_vertex *vtx1, ts_vertex *vtx2, ts_vertex *vtx3){
8f6a69 182     ts_double dX2=vtx2->x-vtx1->x;
SP 183     ts_double dY2=vtx2->y-vtx1->y;
184     ts_double dZ2=vtx2->z-vtx1->z;
185     ts_double dX3=vtx3->x-vtx1->x;
186     ts_double dY3=vtx3->y-vtx1->y;
187     ts_double dZ3=vtx3->z-vtx1->z;
188     ts_double direct=vtx1->x*(dY2*dZ3 -dZ2*dY3)+ 
189         vtx1->y*(dZ2*dX3-dX2*dZ3)+
190         vtx1->z*(dX2*dY3-dY2*dX3);
7958e9 191     return(direct);    
SP 192 }
193
194
195 inline ts_bool vertex_add_tristar(ts_vertex *vtx, ts_triangle *tristarmem){
8f6a69 196     vtx->tristar_no++;
SP 197     vtx->tristar=(ts_triangle **)realloc(vtx->tristar,vtx->tristar_no*sizeof(ts_triangle *));
198     if(vtx->tristar==NULL){
7958e9 199             fatal("Reallocation of memory while adding tristar failed.",3);
SP 200     }
8f6a69 201     vtx->tristar[vtx->tristar_no-1]=tristarmem;
7958e9 202     return TS_SUCCESS;
SP 203 }
204
f74313 205
30ee9c 206
SP 207 /* vtx remove tristar is required in  bondflip. */
208 /* TODO: Check whether it is important to keep the numbering of tristar
209  * elements in some order or not! */
210 inline ts_bool vtx_remove_tristar(ts_vertex *vtx, ts_triangle *tristar){
211     ts_uint i,j=0;
212     for(i=0;i<vtx->tristar_no;i++){
213         if(vtx->tristar[i]!=tristar){
214             vtx->tristar[j]=vtx->tristar[i];
215             j++;
216         }
217     }
218     vtx->tristar_no--;
219     vtx->tristar=realloc(vtx->tristar,vtx->tristar_no*sizeof(ts_triangle *));
220     if(vtx->neigh == NULL){
221             fatal("Reallocation of memory failed during insertion of vertex neighbour in vertex_add_neighbour",3);
222         }
223     return TS_SUCCESS;
224 }
225
226
227
f74313 228 /* ****************************************************************** */
SP 229 /* ***** New vertex copy operations. Inherently they are slow.  ***** */
230 /* ****************************************************************** */
231
b01cc1 232 ts_bool vtx_copy(ts_vertex *cvtx, ts_vertex *ovtx){
f74313 233     memcpy((void *)cvtx,(void *)ovtx,sizeof(ts_vertex));
8ed8fd 234     cvtx->neigh=init_vertex_list(0);
8f6a69 235     cvtx->tristar_no=0;
SP 236     cvtx->bond_no=0;
237     cvtx->tristar=NULL;
238     cvtx->bond=NULL;
239     cvtx->cell=NULL;
b01cc1 240     return TS_SUCCESS;
f74313 241 }
dac2e5 242
SP 243 ts_bool vtx_duplicate(ts_vertex *cvtx, ts_vertex *ovtx){
244     memcpy((void *)cvtx,(void *)ovtx,sizeof(ts_vertex));
245     return TS_SUCCESS;
246 }
247
f74313 248 //TODO: needs to be done
SP 249 ts_vertex **vtx_neigh_copy(ts_vertex_list *vlist,ts_vertex *ovtx){
250         return NULL;
251 }
252
dac2e5 253
SP 254
f74313 255 ts_vertex_list *vertex_list_copy(ts_vertex_list *ovlist){
SP 256     ts_uint i;
257     ts_vertex_list *vlist=(ts_vertex_list *)malloc(sizeof(ts_vertex_list));
b01cc1 258     vlist=memcpy((void *)vlist, (void *)ovlist, sizeof(ts_vertex_list));
f74313 259     ts_vertex **vtx=(ts_vertex **)malloc(vlist->n*sizeof(ts_vertex *));
b01cc1 260     vlist->vtx=vtx;
8f6a69 261     if(vlist->vtx==NULL)
b01cc1 262         fatal("Fatal error reserving memory space for vertex list! Could number of requsted vertices be too large?", 100);
SP 263     for(i=0;i<vlist->n;i++) {
8f6a69 264         vlist->vtx[i]=(ts_vertex *)calloc(1,sizeof(ts_vertex));
b01cc1 265         vlist->vtx[i]->idx=i;
SP 266         vtx_copy(vlist->vtx[i],ovlist->vtx[i]);
f74313 267     }
b01cc1 268
f74313 269     return vlist;
SP 270 }