Trisurf Monte Carlo simulator
Samo Penic
2016-12-06 0634b97e436e480b203a22e1113c473d3ae4d281
commit | author | age
7f6076 1 /* vim: set ts=4 sts=4 sw=4 noet : */
f74313 2 #include "general.h"
d7639a 3 #include<stdio.h>
SP 4 #include "io.h"
f74313 5 #include "vertex.h"
SP 6 #include "bond.h"
d7639a 7 #include<string.h>
a6b1b5 8 #include<stdlib.h>
d7639a 9 #include <sys/types.h>
SP 10 #include <dirent.h>
b14a8d 11 #include "initial_distribution.h"
a2db52 12 #include "poly.h"
055dd7 13 #include "cell.h"
083e03 14 #include <getopt.h>
SP 15 #include <sys/stat.h>
16 #include <sys/types.h>
17 #include <dirent.h>
18 #include <errno.h>
854cb6 19 #include <snapshot.h>
e9eab4 20 /** DUMP STATE TO DISK DRIVE **/
SP 21
1ab449 22 ts_bool dump_state(ts_vesicle *vesicle, ts_uint iteration){
e9eab4 23
SP 24     /* save current state with wrong pointers. Will fix that later */
25     ts_uint i,j,k;
267db5 26     FILE *fh=fopen(command_line_args.dump_fullfilename,"wb");
e9eab4 27
SP 28     /* dump vesicle */
d0a1d3 29     fwrite(vesicle, sizeof(ts_vesicle)-sizeof(ts_double),1,fh);
e9eab4 30     /* dump vertex list */
SP 31     fwrite(vesicle->vlist, sizeof(ts_vertex_list),1,fh);
32     /* dump bond list */
33     fwrite(vesicle->blist, sizeof(ts_bond_list),1,fh);
34     /* dump triangle list */
35     fwrite(vesicle->tlist, sizeof(ts_triangle_list),1,fh);
36     /* dump cell list */
37     fwrite(vesicle->clist, sizeof(ts_cell_list),1,fh);
38     /* dump poly list */
39     fwrite(vesicle->poly_list, sizeof(ts_poly_list),1,fh);
bcf455 40     /* dump filament list */
M 41     fwrite(vesicle->filament_list, sizeof(ts_poly_list),1,fh);
e9eab4 42     /* level 1 complete */
SP 43
44     /*dump vertices*/
45     for(i=0;i<vesicle->vlist->n;i++){
46         fwrite(vesicle->vlist->vtx[i],sizeof(ts_vertex),1,fh);
47         /* dump pointer offsets for:
48                     neigh
49                     bond
50                     tria
51                     cell is ignored
52         */
53         for(j=0;j<vesicle->vlist->vtx[i]->neigh_no;j++){
3c1ac1 54             fwrite(&vesicle->vlist->vtx[i]->neigh[j]->idx,sizeof(ts_uint),1,fh); 
e9eab4 55         }
SP 56         for(j=0;j<vesicle->vlist->vtx[i]->bond_no;j++){
3c1ac1 57             fwrite(&vesicle->vlist->vtx[i]->bond[j]->idx,sizeof(ts_uint),1,fh); 
e9eab4 58         }
SP 59         for(j=0;j<vesicle->vlist->vtx[i]->tristar_no;j++){
3c1ac1 60             fwrite(&vesicle->vlist->vtx[i]->tristar[j]->idx,sizeof(ts_uint),1,fh); 
e9eab4 61         }
SP 62     }
63
64     /*dump bonds*/
65     for(i=0;i<vesicle->blist->n;i++){
66         fwrite(vesicle->blist->bond[i],sizeof(ts_bond),1,fh);
67         /* dump pointer offsets for vtx1 and vtx2 */
3c1ac1 68         //off=(ts_ulong)(vesicle->blist->bond[i]->vtx1-vesicle->vlist->vtx[0]);
SP 69         fwrite(&vesicle->blist->bond[i]->vtx1->idx,sizeof(ts_uint),1,fh); 
70         //off=(ts_ulong)(vesicle->blist->bond[i]->vtx2-vesicle->vlist->vtx[0]);
71         fwrite(&vesicle->blist->bond[i]->vtx2->idx,sizeof(ts_uint),1,fh); 
e9eab4 72     }
SP 73
74     /*dump triangles*/
75     for(i=0;i<vesicle->tlist->n;i++){
76         fwrite(vesicle->tlist->tria[i],sizeof(ts_triangle),1,fh);
77         /* dump pointer offsets for vertex */
3c1ac1 78         fwrite(&vesicle->tlist->tria[i]->vertex[0]->idx,sizeof(ts_uint),1,fh); 
SP 79         fwrite(&vesicle->tlist->tria[i]->vertex[1]->idx,sizeof(ts_uint),1,fh); 
80         fwrite(&vesicle->tlist->tria[i]->vertex[2]->idx,sizeof(ts_uint),1,fh); 
e9eab4 81         /* dump pointer offsets for neigh */
SP 82         for(j=0;j<vesicle->tlist->tria[i]->neigh_no;j++){
3c1ac1 83             fwrite(&vesicle->tlist->tria[i]->neigh[j]->idx,sizeof(ts_uint),1,fh); 
e9eab4 84         }
SP 85     }
86
87
88     /*dump polymeres */
89     for(i=0;i<vesicle->poly_list->n;i++){
90         fwrite(vesicle->poly_list->poly[i],sizeof(ts_poly),1,fh);
3c1ac1 91         fwrite(vesicle->poly_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh);
SP 92         fwrite(vesicle->poly_list->poly[i]->blist,sizeof(ts_bond_list),1,fh);
e9eab4 93     } 
SP 94      
95     /* dump poly vertex(monomer) list*/
96     for(i=0;i<vesicle->poly_list->n;i++){
97         for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
98             fwrite(vesicle->poly_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh);
99             /* dump offset for neigh and bond */
100             for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh_no;k++){
3c1ac1 101                // off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]-vesicle->poly_list->poly[i]->vlist->vtx[0]);
SP 102                 fwrite(&vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]->idx,sizeof(ts_uint),1,fh); 
e9eab4 103             }
SP 104             for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->bond_no;k++){
3c1ac1 105                 //off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]-vesicle->poly_list->poly[i]->blist->bond[0]);
SP 106                 fwrite(&vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]->idx,sizeof(ts_uint),1,fh); 
e9eab4 107             }
SP 108         }
3c1ac1 109     // grafted vtx on vesicle data dump
SP 110         fwrite(&vesicle->poly_list->poly[i]->grafted_vtx->idx, sizeof(ts_uint),1,fh);
e9eab4 111     }
SP 112     /* dump poly bonds between monomers list*/
113     for(i=0;i<vesicle->poly_list->n;i++){
114         for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
115             fwrite(vesicle->poly_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh);
116             /* dump vtx1 and vtx2 offsets */
3c1ac1 117             //off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx1-vesicle->poly_list->poly[i]->vlist->vtx[0]);
SP 118             fwrite(&vesicle->poly_list->poly[i]->blist->bond[j]->vtx1->idx,sizeof(ts_uint),1,fh); 
119 //            off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx2-vesicle->poly_list->poly[i]->vlist->vtx[0]);
120             fwrite(&vesicle->poly_list->poly[i]->blist->bond[j]->vtx2->idx,sizeof(ts_uint),1,fh); 
e9eab4 121         }
SP 122     }
123
bcf455 124
M 125   /*dump filamentes grandes svinjas */
126     for(i=0;i<vesicle->filament_list->n;i++){
127         fwrite(vesicle->filament_list->poly[i],sizeof(ts_poly),1,fh);
128         fwrite(vesicle->filament_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh);
129         fwrite(vesicle->filament_list->poly[i]->blist,sizeof(ts_bond_list),1,fh);
130     } 
131      
132     /* dump filamentes vertex(monomer) list*/
133     for(i=0;i<vesicle->filament_list->n;i++){
134         for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){
135             fwrite(vesicle->filament_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh);
136             /* dump offset for neigh and bond */
137             for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh_no;k++){
138                // off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]-vesicle->poly_list->poly[i]->vlist->vtx[0]);
139                 fwrite(&vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh[k]->idx,sizeof(ts_uint),1,fh); 
140             }
141             for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->bond_no;k++){
142                 //off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]-vesicle->poly_list->poly[i]->blist->bond[0]);
143                 fwrite(&vesicle->filament_list->poly[i]->vlist->vtx[j]->bond[k]->idx,sizeof(ts_uint),1,fh); 
144             }
145         }
146     }
147     /* dump poly bonds between monomers list*/
148     for(i=0;i<vesicle->filament_list->n;i++){
149         for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){
150             fwrite(vesicle->filament_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh);
151             /* dump vtx1 and vtx2 offsets */
152             //off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx1-vesicle->poly_list->poly[i]->vlist->vtx[0]);
153             fwrite(&vesicle->filament_list->poly[i]->blist->bond[j]->vtx1->idx,sizeof(ts_uint),1,fh); 
154 //            off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx2-vesicle->poly_list->poly[i]->vlist->vtx[0]);
155             fwrite(&vesicle->filament_list->poly[i]->blist->bond[j]->vtx2->idx,sizeof(ts_uint),1,fh); 
156         }
157     }
158
159
160
e9eab4 161 /* pointer offsets for fixing the restored pointers */
SP 162 /* need pointers for 
163     vlist->vtx
164     blist->bond
165     tlist->tria
166     clist->cell
167     poly_list->poly
168     and for each poly:
169         poly_list->poly->vtx
170         poly_list->poly->bond
171 */
172
bd2210 173 //    fwrite(vesicle->clist, sizeof(ts_cell_list),1,  fh);
9166cb 174 /* write tape information on vesicle */
bd2210 175 //    fwrite(vesicle->tape,sizeof(ts_tape),1,fh);
1ab449 176     fwrite(&iteration, sizeof(ts_uint),1,fh);
e9eab4 177     fclose(fh);
SP 178     return TS_SUCCESS;
179 }
180
181
182 /** RESTORE DUMP FROM DISK **/
1ab449 183 ts_vesicle *restore_state(ts_uint *iteration){
e9eab4 184     ts_uint i,j,k;
267db5 185     FILE *fh=fopen(command_line_args.dump_fullfilename,"rb");
3f5c83 186
SP 187     struct stat sb;
267db5 188     if (stat(command_line_args.dump_fullfilename, &sb) == -1) {
3f5c83 189         //dump file does not exist.
SP 190         return NULL;
191     }
192
193     //check if it is regular file
194     if((sb.st_mode & S_IFMT) != S_IFREG) {
195         //dump file is not a regular file.
196         ts_fprintf(stderr,"Dump file is not a regular file!\n");
197         return NULL;
198     }
199
e9eab4 200     ts_uint retval;
3c1ac1 201     ts_uint idx;
e9eab4 202
SP 203 /* we restore all the data from the dump */
204     /* restore vesicle */
205     ts_vesicle *vesicle=(ts_vesicle *)calloc(1,sizeof(ts_vesicle));
d0a1d3 206     retval=fread(vesicle, sizeof(ts_vesicle)-sizeof(ts_double),1,fh);
86f5e7 207 //    fprintf(stderr,"was here! %e\n",vesicle->dmax);
SP 208
e9eab4 209     /* restore vertex list */
SP 210     vesicle->vlist=(ts_vertex_list *)malloc(sizeof(ts_vertex_list));
211     retval=fread(vesicle->vlist, sizeof(ts_vertex_list),1,fh);
212     /* restore bond list */
213     vesicle->blist=(ts_bond_list *)malloc(sizeof(ts_bond_list));
214     retval=fread(vesicle->blist, sizeof(ts_bond_list),1,fh);
215     /* restore triangle list */
216     vesicle->tlist=(ts_triangle_list *)malloc(sizeof(ts_triangle_list));
217     retval=fread(vesicle->tlist, sizeof(ts_triangle_list),1,fh);
218     /* restore cell list */
219     vesicle->clist=(ts_cell_list *)malloc(sizeof(ts_cell_list));
220     retval=fread(vesicle->clist, sizeof(ts_cell_list),1,fh);
221     /* restore poly list */
222     vesicle->poly_list=(ts_poly_list *)calloc(1,sizeof(ts_poly_list));
223     retval=fread(vesicle->poly_list, sizeof(ts_poly_list),1,fh);
bcf455 224     /* restore filament list */
M 225     vesicle->filament_list=(ts_poly_list *)calloc(1,sizeof(ts_poly_list));
226     retval=fread(vesicle->filament_list, sizeof(ts_poly_list),1,fh);
e9eab4 227     /* level 1 complete */
SP 228
3c1ac1 229 /* prerequisity. Bonds must be malloced before vertexes are recreated */
SP 230   vesicle->blist->bond=(ts_bond **)calloc(vesicle->blist->n,sizeof(ts_bond *));
231     for(i=0;i<vesicle->blist->n;i++){
232         vesicle->blist->bond[i]=(ts_bond *)malloc(sizeof(ts_bond));
233     }
234 /* prerequisity. Triangles must be malloced before vertexes are recreated */
235   vesicle->tlist->tria=(ts_triangle **)calloc(vesicle->tlist->n,sizeof(ts_triangle *));
236     for(i=0;i<vesicle->tlist->n;i++){
237         vesicle->tlist->tria[i]=(ts_triangle *)malloc(sizeof(ts_triangle));
238 }
a37ba2 239 /* prerequisity. Vertices must be malloced before vertexes are recreated */
e9eab4 240     vesicle->vlist->vtx=(ts_vertex **)calloc(vesicle->vlist->n,sizeof(ts_vertex *));
a37ba2 241  for(i=0;i<vesicle->vlist->n;i++){
e9eab4 242         vesicle->vlist->vtx[i]=(ts_vertex *)malloc(sizeof(ts_vertex));
a37ba2 243  }
SP 244   /*restore vertices*/
245     for(i=0;i<vesicle->vlist->n;i++){
e9eab4 246         retval=fread(vesicle->vlist->vtx[i],sizeof(ts_vertex),1,fh);
SP 247         /*restore neigh, bond, tristar. Ignoring cell */
248         vesicle->vlist->vtx[i]->neigh=(ts_vertex **)calloc(vesicle->vlist->vtx[i]->neigh_no, sizeof(ts_vertex *));
249         for(j=0;j<vesicle->vlist->vtx[i]->neigh_no;j++){
3c1ac1 250             retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 251             vesicle->vlist->vtx[i]->neigh[j]=vesicle->vlist->vtx[idx];
e9eab4 252         }
SP 253         vesicle->vlist->vtx[i]->bond=(ts_bond **)calloc(vesicle->vlist->vtx[i]->bond_no, sizeof(ts_bond *));
254         for(j=0;j<vesicle->vlist->vtx[i]->bond_no;j++){
3c1ac1 255             retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 256 /* pointer can be assigned only when list of bonds is fully initialized in memory. Thus bondlist popularization must be done before vertex can reference to it */
257             vesicle->vlist->vtx[i]->bond[j]=vesicle->blist->bond[idx];    
e9eab4 258         }
SP 259
260         vesicle->vlist->vtx[i]->tristar=(ts_triangle **)calloc(vesicle->vlist->vtx[i]->tristar_no, sizeof(ts_triangle *));
261         for(j=0;j<vesicle->vlist->vtx[i]->tristar_no;j++){
3c1ac1 262             retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 263 /* same comment as above */
264             vesicle->vlist->vtx[i]->tristar[j]=vesicle->tlist->tria[idx];
e9eab4 265         }
SP 266
267     }
268
269     /*restore bonds*/
3c1ac1 270    // vesicle->blist->bond=(ts_bond **)calloc(vesicle->blist->n,sizeof(ts_bond *)); // done before.
e9eab4 271     for(i=0;i<vesicle->blist->n;i++){
3c1ac1 272      //   vesicle->blist->bond[i]=(ts_bond *)malloc(sizeof(ts_bond)); //done before.
e9eab4 273         retval=fread(vesicle->blist->bond[i],sizeof(ts_bond),1,fh);
SP 274         /* restore vtx1 and vtx2 */
3c1ac1 275         retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 276         vesicle->blist->bond[i]->vtx1=vesicle->vlist->vtx[idx];
277         retval=fread(&idx,sizeof(ts_uint),1,fh);
278         vesicle->blist->bond[i]->vtx2=vesicle->vlist->vtx[idx];
e9eab4 279     }
SP 280
281     /*restore triangles*/
3c1ac1 282 //    vesicle->tlist->tria=(ts_triangle **)calloc(vesicle->tlist->n,sizeof(ts_triangle *)); // done before
e9eab4 283     for(i=0;i<vesicle->tlist->n;i++){
3c1ac1 284  //       vesicle->tlist->tria[i]=(ts_triangle *)malloc(sizeof(ts_triangle)); // done before
e9eab4 285         retval=fread(vesicle->tlist->tria[i],sizeof(ts_triangle),1,fh);
SP 286         /* restore pointers for vertices */
3c1ac1 287         retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 288         vesicle->tlist->tria[i]->vertex[0]=vesicle->vlist->vtx[idx];
289         retval=fread(&idx,sizeof(ts_uint),1,fh);
290         vesicle->tlist->tria[i]->vertex[1]=vesicle->vlist->vtx[idx];
291         retval=fread(&idx,sizeof(ts_uint),1,fh);
292         vesicle->tlist->tria[i]->vertex[2]=vesicle->vlist->vtx[idx];
e9eab4 293         /* restore pointers for neigh */
3c1ac1 294      vesicle->tlist->tria[i]->neigh=(ts_triangle **)malloc(vesicle->tlist->tria[i]->neigh_no*sizeof(ts_triangle *));
e9eab4 295         for(j=0;j<vesicle->tlist->tria[i]->neigh_no;j++){
3c1ac1 296             retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 297             vesicle->tlist->tria[i]->neigh[j]=vesicle->tlist->tria[idx];
e9eab4 298         }
SP 299
300     }
301    
302     /*restore cells */
303 /*TODO: do we need to recalculate cells here? */
304 /*    vesicle->clist->cell=(ts_cell **)malloc(vesicle->clist->cellno*sizeof(ts_cell *));
305     for(i=0;i<vesicle->clist->cellno;i++){
306         vesicle->clist->cell[i]=(ts_cell *)malloc(sizeof(ts_cell));
307         retval=fread(vesicle->clist->cell[i],sizeof(ts_cell),1,fh);
308     }
309 */
310     /*restore polymeres */
311     vesicle->poly_list->poly = (ts_poly **)calloc(vesicle->poly_list->n,sizeof(ts_poly *));
312     for(i=0;i<vesicle->poly_list->n;i++){
313         vesicle->poly_list->poly[i]=(ts_poly *)calloc(1,sizeof(ts_poly));
314         retval=fread(vesicle->poly_list->poly[i],sizeof(ts_poly),1,fh);
3c1ac1 315         vesicle->poly_list->poly[i]->vlist=(ts_vertex_list *)calloc(1,sizeof(ts_vertex_list));
SP 316         retval=fread(vesicle->poly_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh);
317         vesicle->poly_list->poly[i]->blist=(ts_bond_list *)calloc(1,sizeof(ts_bond_list));
318         retval=fread(vesicle->poly_list->poly[i]->blist,sizeof(ts_bond_list),1,fh);
319     /* initialize adress space for pointers that will hold specific vertices (monomers) and bonds */
320         vesicle->poly_list->poly[i]->vlist->vtx=(ts_vertex **)calloc(vesicle->poly_list->poly[i]->vlist->n,sizeof(ts_vertex *));
321         vesicle->poly_list->poly[i]->blist->bond=(ts_bond **)calloc(vesicle->poly_list->poly[i]->blist->n,sizeof(ts_bond *));
322      for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
323             vesicle->poly_list->poly[i]->vlist->vtx[j]=(ts_vertex *)malloc(sizeof(ts_vertex));
324     }
325     for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
326             vesicle->poly_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond));
327     }
328
e9eab4 329     } 
3c1ac1 330
e9eab4 331      
SP 332     /* restore poly vertex(monomer) list*/
333     for(i=0;i<vesicle->poly_list->n;i++){
334         for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
335             retval=fread(vesicle->poly_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh);
3c1ac1 336                 
e9eab4 337             /* restore neigh and bonds */
SP 338             vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh=(ts_vertex **)calloc(vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh_no, sizeof(ts_vertex *));
339             for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh_no;k++){
3c1ac1 340                 retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 341                 vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]=vesicle->poly_list->poly[i]->vlist->vtx[idx];
e9eab4 342             }
SP 343             vesicle->poly_list->poly[i]->vlist->vtx[j]->bond=(ts_bond **)calloc(vesicle->poly_list->poly[i]->vlist->vtx[j]->bond_no, sizeof(ts_bond *));
344             for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->bond_no;k++){
3c1ac1 345                 retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 346                 vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]=vesicle->poly_list->poly[i]->blist->bond[idx];
e9eab4 347             }
SP 348
349         }
3c1ac1 350     /* restore grafted vtx on vesicle and grafted_poly */
SP 351                 retval=fread(&idx,sizeof(ts_uint),1,fh);
352         vesicle->vlist->vtx[idx]->grafted_poly=vesicle->poly_list->poly[i];
353         vesicle->poly_list->poly[i]->grafted_vtx=vesicle->vlist->vtx[idx];    
e9eab4 354     }
3c1ac1 355
e9eab4 356     /* restore poly bonds between monomers list*/
SP 357     for(i=0;i<vesicle->poly_list->n;i++){
358         for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
3c1ac1 359        //     vesicle->poly_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond));
e9eab4 360             retval=fread(vesicle->poly_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh);
SP 361             /* restore vtx1 and vtx2 */
3c1ac1 362                 retval=fread(&idx,sizeof(ts_uint),1,fh);
SP 363                 vesicle->poly_list->poly[i]->blist->bond[j]->vtx1=vesicle->poly_list->poly[i]->vlist->vtx[idx];
364                 retval=fread(&idx,sizeof(ts_uint),1,fh);
365                 vesicle->poly_list->poly[i]->blist->bond[j]->vtx2=vesicle->poly_list->poly[i]->vlist->vtx[idx];
e9eab4 366         }
SP 367     }
368
bcf455 369     /*restore filaments */
M 370     vesicle->filament_list->poly = (ts_poly **)calloc(vesicle->filament_list->n,sizeof(ts_poly *));
371     for(i=0;i<vesicle->filament_list->n;i++){
372         vesicle->filament_list->poly[i]=(ts_poly *)calloc(1,sizeof(ts_poly));
373         retval=fread(vesicle->filament_list->poly[i],sizeof(ts_poly),1,fh);
374         vesicle->filament_list->poly[i]->vlist=(ts_vertex_list *)calloc(1,sizeof(ts_vertex_list));
375         retval=fread(vesicle->filament_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh);
376         vesicle->filament_list->poly[i]->blist=(ts_bond_list *)calloc(1,sizeof(ts_bond_list));
377         retval=fread(vesicle->filament_list->poly[i]->blist,sizeof(ts_bond_list),1,fh);
378     /* initialize adress space for pointers that will hold specific vertices (monomers) and bonds */
379         vesicle->filament_list->poly[i]->vlist->vtx=(ts_vertex **)calloc(vesicle->filament_list->poly[i]->vlist->n,sizeof(ts_vertex *));
380         vesicle->filament_list->poly[i]->blist->bond=(ts_bond **)calloc(vesicle->filament_list->poly[i]->blist->n,sizeof(ts_bond *));
381      for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){
382             vesicle->filament_list->poly[i]->vlist->vtx[j]=(ts_vertex *)malloc(sizeof(ts_vertex));
383     }
384     for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){
385             vesicle->filament_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond));
386     }
387
388     } 
389
390      
391     /* restore poly vertex(monomer) list*/
392     for(i=0;i<vesicle->filament_list->n;i++){
393         for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){
394             retval=fread(vesicle->filament_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh);
395                 
396             /* restore neigh and bonds */
397             vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh=(ts_vertex **)calloc(vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh_no, sizeof(ts_vertex *));
398             for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh_no;k++){
399                 retval=fread(&idx,sizeof(ts_uint),1,fh);
400                 vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh[k]=vesicle->filament_list->poly[i]->vlist->vtx[idx];
401             }
402             vesicle->filament_list->poly[i]->vlist->vtx[j]->bond=(ts_bond **)calloc(vesicle->filament_list->poly[i]->vlist->vtx[j]->bond_no, sizeof(ts_bond *));
403             for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->bond_no;k++){
404                 retval=fread(&idx,sizeof(ts_uint),1,fh);
405                 vesicle->filament_list->poly[i]->vlist->vtx[j]->bond[k]=vesicle->filament_list->poly[i]->blist->bond[idx];
406             }
407
408         }
409     }
410
411     /* restore poly bonds between monomers list*/
412     for(i=0;i<vesicle->filament_list->n;i++){
413         for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){
414        //     vesicle->poly_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond));
415             retval=fread(vesicle->filament_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh);
416             /* restore vtx1 and vtx2 */
417                 retval=fread(&idx,sizeof(ts_uint),1,fh);
418                 vesicle->filament_list->poly[i]->blist->bond[j]->vtx1=vesicle->filament_list->poly[i]->vlist->vtx[idx];
419                 retval=fread(&idx,sizeof(ts_uint),1,fh);
420                 vesicle->filament_list->poly[i]->blist->bond[j]->vtx2=vesicle->filament_list->poly[i]->vlist->vtx[idx];
421         }
422     }
267db5 423     vesicle->tape=parsetape(command_line_args.tape_fullfilename);
063289 424 // recreating space for cells // 
055dd7 425     vesicle->clist=init_cell_list(vesicle->tape->ncxmax, vesicle->tape->ncymax, vesicle->tape->nczmax, vesicle->tape->stepsize);
364c30 426     vesicle->clist->max_occupancy=8;
47a7ac 427 //    vesicle->tape=(ts_tape *)malloc(sizeof(ts_tape));
SP 428 //    retval=fread(vesicle->tape, sizeof(ts_tape),1,fh);
1ab449 429     retval=fread(iteration,sizeof(ts_uint),1,fh);
e9eab4 430     if(retval); 
SP 431     fclose(fh);
432     return vesicle;
433 }
434
435
436
083e03 437 ts_bool parse_args(int argc, char **argv){
3f5c83 438     int c, retval;
SP 439     struct stat sb;
440     sprintf(command_line_args.path, "./"); //clear string;
7b0c07 441     sprintf(command_line_args.output_fullfilename,"output.pvd");
SP 442     sprintf(command_line_args.dump_fullfilename,"dump.bin");
443     sprintf(command_line_args.tape_fullfilename,"tape");
1bf3c3 444     sprintf(command_line_args.tape_templatefull,"./tape");
3f5c83 445             FILE *file;
SP 446     
083e03 447 while (1)
SP 448      {
449        static struct option long_options[] =
450          {
8a6614 451            {"force-from-tape", no_argument,       &(command_line_args.force_from_tape), 1},
SP 452        {"reset-iteration-count", no_argument, &(command_line_args.reset_iteration_count), 1},
083e03 453            {"tape",     no_argument,       0, 't'},
fd8126 454        {"version", no_argument, 0, 'v'},
083e03 455            {"output-file",  required_argument, 0, 'o'},
SP 456            {"directory",  required_argument, 0, 'd'},
3f5c83 457            {"dump-filename", required_argument,0, 'f'},
7b0c07 458            {"tape-options",required_argument,0,'c'},
1bf3c3 459            {"tape-template", required_argument,0,0},
ab798b 460             {"restore-from-vtk",required_argument,0,0},
083e03 461            {0, 0, 0, 0}
SP 462          };
463        /* getopt_long stores the option index here. */
464        int option_index = 0;
465
fd8126 466        c = getopt_long (argc, argv, "d:f:o:t:c:v",
083e03 467                         long_options, &option_index);
SP 468
469        /* Detect the end of the options. */
470        if (c == -1)
471          break;
472
473        switch (c)
474          {
475          case 0:
476            /* If this option set a flag, do nothing else now. */
477            if (long_options[option_index].flag != 0)
478              break;
1bf3c3 479 /*           printf ("option %s", long_options[option_index].name);
083e03 480            if (optarg)
1bf3c3 481              printf (" with arg %s", optarg); 
SP 482            printf ("\n"); */
483             //TODO: find a better way.
484             if(strcmp(long_options[option_index].name,"tape-template")==0){
485                 strcpy(command_line_args.tape_templatefull,optarg);
486             }
ab798b 487             if(strcmp(long_options[option_index].name,"restore-from-vtk")==0){
ee84bd 488                 strcpy(command_line_args.dump_from_vtk,optarg);
SP 489             }
083e03 490            break;
fd8126 491      case 'v':
SP 492         fprintf(stdout,"TRISURF-NG v. %s, compiled on: %s %s.\n", TS_VERSION, __DATE__, __TIME__);
493             fprintf(stdout,"Programming done by: Samo Penic and Miha Fosnaric\n");
494             fprintf(stdout,"Released under terms of GPLv3\n");
495         exit(0);
083e03 496
7b0c07 497          case 'c':
SP 498               strcpy(command_line_args.tape_opts,optarg);
499             break;
500          case 't': //tape
3f5c83 501                 strcpy(command_line_args.tape_fullfilename,optarg);
083e03 502            break;
SP 503
7b0c07 504          case 'o':  //set filename of master pvd output file
3f5c83 505             strcpy(command_line_args.output_fullfilename, optarg);
SP 506             break;
083e03 507
SP 508          case 'd':
509             //check if directory exists. If not create one. If creation is
510             //successful, set directory for output files.
511             //printf ("option -d with value `%s'\n", optarg);
3f5c83 512             if (stat(optarg, &sb) == -1) {
SP 513                 //directory does not exist
083e03 514                 retval=mkdir(optarg, 0700);
SP 515                 if(retval){
3f5c83 516                     fatal("Could not create requested directory. Check if you have permissions",1);
083e03 517                 }
SP 518             }
3f5c83 519             //check if is a proper directory
SP 520             else if((sb.st_mode & S_IFMT) != S_IFDIR) {
521                 //it is not a directory. fatal error.
522                 ts_fprintf(stderr,"%s is not a directory!\n",optarg);
523                 fatal("Cannot continue",1);
083e03 524             }
3f5c83 525             strcpy(command_line_args.path, optarg);
083e03 526            break;
SP 527
528         case 'f':
3f5c83 529             strcpy(command_line_args.dump_fullfilename, optarg);
083e03 530             break;
SP 531
532          case '?':
533            /* getopt_long already printed an error message. */
ba73ab 534             print_help(stdout);
SP 535 //ts_fprintf(stderr,"\n\nhere comes the help.\n\n");
083e03 536             fatal("Ooops, read help first",1);
SP 537            break;
538
539          default:
540            exit (1);
541          }
542      }
543
7b0c07 544 //Here we set correct values for full filenames!
SP 545     char *buffer=(char *)malloc(10000*sizeof(char));
546     //correct the path and add trailing /
547     if(command_line_args.path[strlen(command_line_args.path)-1]!='/') strcat(command_line_args.path,"/");
548    
549 /* master pvd output file */ 
550     strcpy(buffer,command_line_args.path);
551     strcat(buffer,command_line_args.output_fullfilename);
552     if ((file = fopen(buffer, "w")) == NULL) {
553                 fprintf(stderr,"Could not create output file %s!\n", buffer);
554                 fatal("Please specify correct output file or check permissions of the file",1);
1bf3c3 555                 //there is a tape template. make a copy into desired directory
7b0c07 556                 
SP 557             } else {
558                 fclose(file);
559                 strcpy(command_line_args.output_fullfilename,buffer);
560             }
561
562 /* tape file */
563     strcpy(buffer,command_line_args.path);
564     strcat(buffer,command_line_args.tape_fullfilename);
565     if (stat(buffer, &sb) == -1) {
1bf3c3 566
SP 567                 //tape does not exist. does tape template exist?
568                 if(stat(command_line_args.tape_templatefull, &sb)==-1){ 
569                     ts_fprintf(stderr,"Tape '%s' does not exist and no tape template was specified (or does not exist)!\n",buffer);
570                     fatal("Please select correct tape or check permissions of the file",1);
571                 } else {
572                     //tape template found
573                     fatal("Samo did not program template copy yet",1); 
574                 }
7b0c07 575             } else {
SP 576                 strcpy(command_line_args.tape_fullfilename,buffer);
577             }
578
579
580 /* dump file */
581             strcpy(buffer,command_line_args.path);
582             strcat(buffer,command_line_args.dump_fullfilename);
583             //check if dump file exist first.
584             if (stat(buffer, &sb) == -1) {
585                 //no dump file. check if we can create one.
586                 if ((file = fopen(buffer, "w")) == NULL) {
587                     fprintf(stderr,"Could not create dump file '%s'!\n",buffer);
588                     fatal("Please specify correct dump file or check permissions of the file",1);
589                 } else {
590                     fclose(file);
591                     //good, file is writeable. delete it for now.
592                     remove(buffer);
593                 }
594             }
595             strcpy(command_line_args.dump_fullfilename, buffer);
596
597
598     free(buffer);
083e03 599     return TS_SUCCESS;
SP 600
601 }
602
603
ba73ab 604 ts_bool print_help(FILE *fd){
SP 605     fprintf(fd,"TRISURF-NG v. %s, compiled on: %s %s.\n", TS_VERSION, __DATE__, __TIME__);
606     fprintf(fd,"Programming done by: Samo Penic and Miha Fosnaric\n");
607     fprintf(fd,"Released under terms of GPLv3\n\n");
608
609     fprintf(fd, "Invoking trisurf-ng without any flags results in default operation. Program reads 'tape' file and 'dump.bin' binary representation of the vesicle from disk and continues the simulation where it was aborted (as specified in 'dump.bin').\n\n");
610     fprintf(fd,"If 'tape' has different values than binary dump, those are used (if possible -- some values in tape, such as nshell cannot be modified).\n\n\n");
611     fprintf(fd,"However, if dump.bin is not present, user is notified to specify --force-from-tape flag. The vesicle will be created from specifications in tape.\n\n\n");
612     fprintf(fd,"Flags:\n\n");
613     fprintf(fd,"--force-from-tape\t\t makes initial shape of the vesicle from tape. Ignores already existing binary dump and possible simulation results.\n");
614     fprintf(fd,"--restore-from-vtk\t\t VTK's file ending with '.vtu' are preferred way to make state snapshots for restoration. With this flag the restoration of the vesicle from vtk is possible. The simulation will continue if hidden '.status' file with last iteration done is present. Otherwise it will start simulation from timestep 0.\n");
c4178c 615     fprintf(fd,"--reset-iteration-count\t\t starts simulation from the beginning (using binary dump).\n");
ba73ab 616     fprintf(fd,"--tape (or -t)\t\t specifies tape filename. For --force-from-tape and restoring from binary dump. Defaults to 'tape'.\n");
SP 617     fprintf(fd,"--version (or -v)\t\t Prints version information.\n");
618     fprintf(fd,"--output-file (or -o)\t\t Specifies filename of .PVD file. Defaults to 'output.pvd'\n");
c4178c 619     fprintf(fd,"--dump-filename (or -f)\t\t specifies filename for binary dump&restore. Defaults to 'dump.bin'\n\n\n");
SP 620     fprintf(fd,"Examples:\n\n");
621     fprintf(fd,"trisurf --force-from-tape\n");
622     fprintf(fd,"trisurf --reset-iteration-count\n");
623     fprintf(fd,"trisurf --restore-from-vtk filename.vtu\n");
624     fprintf(fd,"\n\n");
625
ba73ab 626     return TS_SUCCESS;
SP 627 }
628
629
083e03 630
d7639a 631 ts_bool print_vertex_list(ts_vertex_list *vlist){
SP 632     ts_uint i;
633     printf("Number of vertices: %u\n",vlist->n);
634     for(i=0;i<vlist->n;i++){
a6b1b5 635         printf("%u: %f %f %f\n",
8f6a69 636 vlist->vtx[i]->idx,vlist->vtx[i]->x, vlist->vtx[i]->y, vlist->vtx[i]->z);
d7639a 637     }
SP 638     return TS_SUCCESS;
639 }
640
641 ts_bool print_vertex_neighbours(ts_vertex_list *vlist){
642     ts_uint i,j;
a6b1b5 643     ts_vertex **vtx=vlist->vtx;
d7639a 644     printf("Vertex id(neigh no): (neighvertex coord) (neighvertex coord) ...\n");
SP 645     for(i=0;i<vlist->n;i++){
8f6a69 646         printf("%u(%u): ",vtx[i]->idx,vtx[i]->neigh_no);
SP 647         for(j=0;j<vtx[i]->neigh_no;j++){
648             printf("(%f,%f,%f)",vtx[i]->neigh[j]->x,
649 vtx[i]->neigh[j]->y,vtx[i]->neigh[j]->z);
d7639a 650         }
SP 651         printf("\n");
652     }
653
654 return TS_SUCCESS;
655 }
656
657 ts_bool write_vertex_fcompat_file(ts_vertex_list *vlist,ts_char *filename){
a6b1b5 658     ts_vertex **vtx=vlist->vtx;
d7639a 659     ts_uint i;
SP 660     FILE *fh;
661     
662     fh=fopen(filename, "w");
663     if(fh==NULL){
664         err("Cannot open file %s for writing");
665         return TS_FAIL;
666     }
667     for(i=0;i<vlist->n;i++)
8f6a69 668         fprintf(fh," %E\t%E\t%E\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
d7639a 669
SP 670     fclose(fh);
671 return TS_SUCCESS;
672 }
673
674
675 ts_bool fprint_vertex_list(FILE *fh,ts_vertex_list *vlist){
676     ts_uint i,j;
677     for(i=0;i<vlist->n;i++){
8f6a69 678         fprintf(fh," %.17E\t%.17E\t%.17E\t%u\n",vlist->vtx[i]->x,
SP 679             vlist->vtx[i]->y, vlist->vtx[i]->z,
680             vlist->vtx[i]->neigh_no);
681         for(j=0;j<vlist->vtx[i]->neigh_no;j++){
682             fprintf(fh,"\t%u",(ts_uint)(vlist->vtx[i]->neigh[j]->idx));
a6b1b5 683         //-vlist->vtx+1));
d7639a 684         }
SP 685         fprintf(fh,"\n");
686     }
687     return TS_SUCCESS;
688 }
689
690 ts_bool fprint_tristar(FILE *fh, ts_vesicle *vesicle){
691     ts_uint i,j;
a6b1b5 692     for(i=0;i<vesicle->vlist->n;i++){
8f6a69 693         fprintf(fh,"\t%u",vesicle->vlist->vtx[i]->tristar_no);
SP 694         for(j=0;j<vesicle->vlist->vtx[i]->tristar_no;j++){
695             fprintf(fh,"\t%u",(ts_uint)(vesicle->vlist->vtx[i]->tristar[j]->idx));//-vesicle->tlist->tria+1));
d7639a 696         }
SP 697         fprintf(fh,"\n");
698     }
699     return TS_SUCCESS;
700 }
701
702 ts_bool fprint_triangle_list(FILE *fh, ts_vesicle *vesicle){
a6b1b5 703         ts_triangle_list *tlist=vesicle->tlist;
d7639a 704       ts_uint i,j;
SP 705     for(i=0;i<tlist->n;i++){
41a035 706         fprintf(fh,"\t%u",tlist->tria[i]->neigh_no);
SP 707         for(j=0;j<tlist->tria[i]->neigh_no;j++){
708             fprintf(fh,"\t%u",(ts_uint)(tlist->tria[i]->neigh[j]->idx));//-tlist->tria+1)); 
d7639a 709         }
SP 710         fprintf(fh,"\n");
711             for(j=0;j<3;j++){
41a035 712             fprintf(fh,"\t%u",(ts_uint)(tlist->tria[i]->vertex[j]->idx));//-vesicle->vlist->vtx+1)); 
d7639a 713             }
SP 714         fprintf(fh,"\n");
41a035 715         fprintf(fh,"%.17E\t%.17E\t%.17E\n",tlist->tria[i]->xnorm,
SP 716 tlist->tria[i]->ynorm,tlist->tria[i]->znorm);
d7639a 717         fprintf(fh,"0.00000000000000000\n0.00000000000000000\n");
SP 718     }
719     return TS_SUCCESS;
720 }
721
722 ts_bool fprint_vertex_data(FILE *fh,ts_vertex_list *vlist){
723     ts_uint i,j;
724     for(i=0;i<vlist->n;i++){
725         fprintf(fh," %.17E\t%.17E\t%.17E\t%.17E\t%.17E\t%u\n",
8f6a69 726         vlist->vtx[i]->xk,vlist->vtx[i]->c,vlist->vtx[i]->energy,
SP 727         vlist->vtx[i]->energy_h, vlist->vtx[i]->curvature, 0);
cef6ca 728         for(j=0;j<vlist->vtx[i]->bond_no;j++){
8f6a69 729             fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length_dual);
d7639a 730         }
SP 731             fprintf(fh,"\n");
cef6ca 732         for(j=0;j<vlist->vtx[i]->bond_no;j++){
8f6a69 733             fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length);
d7639a 734         }
SP 735             fprintf(fh,"\n");
736     }
737     return TS_SUCCESS;
738 }
739
740 ts_bool fprint_bonds(FILE *fh,ts_vesicle *vesicle){
741     ts_uint i;
a6b1b5 742     for(i=0;i<vesicle->blist->n;i++){
e016c4 743         fprintf(fh,"\t%u\t%u\n",(ts_uint)(vesicle->blist->bond[i]->vtx1->idx),
a6b1b5 744 //-vesicle->vlist->vtx+1),
e016c4 745         (ts_uint)(vesicle->blist->bond[i]->vtx2->idx));
a6b1b5 746     //-vesicle->vlist.vtx+1));
d7639a 747     }
SP 748     return TS_SUCCESS;
749 }
750
751
752 ts_bool write_dout_fcompat_file(ts_vesicle *vesicle, ts_char *filename){
753     FILE *fh;
754     fh=fopen(filename, "w");
755     if(fh==NULL){
756         err("Cannot open file %s for writing");
757         return TS_FAIL;
758     }
759     fprintf(fh,"%.17E\n%.17E\n",vesicle->stepsize,vesicle->dmax);
a6b1b5 760     fprint_vertex_list(fh,vesicle->vlist);
d7639a 761     fprint_tristar(fh,vesicle);
SP 762     fprint_triangle_list(fh,vesicle);
a6b1b5 763     fprint_vertex_data(fh,vesicle->vlist);
d7639a 764     fprint_bonds(fh,vesicle);
SP 765     fclose(fh);    
766     return TS_SUCCESS;
767 }
768
769 ts_bool read_tape_fcompat_file(ts_vesicle *vesicle, ts_char *filename){
770     FILE *fh;
771     char line[255];
772     fh=fopen(filename, "r");
773         if(fh==NULL){
774                 err("Cannot open file for reading... Nonexistant file?");
775                 return TS_FAIL;
776         }
a6b1b5 777     ts_uint retval=1;
d7639a 778     while(retval!=EOF){
SP 779         retval=fscanf(fh,"%s",line);
780         
781         fprintf(stderr,"%s",line);
782     }    
783     fclose(fh);    
784     return TS_SUCCESS;
785 }
786
787 ts_bool write_master_xml_file(ts_char *filename){
788      FILE *fh;
789     ts_char *i,*j;
790     ts_uint tstep;
791     ts_char *number;
792         fh=fopen(filename, "w");
793         if(fh==NULL){
794                 err("Cannot open file %s for writing");
795                 return TS_FAIL;
796         }
797
798     fprintf(fh,"<?xml version=\"1.0\"?>\n<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n<Collection>");
267db5 799     DIR *dir = opendir(command_line_args.path);
d7639a 800     if(dir){
SP 801         struct dirent *ent;
802         tstep=0;
803         while((ent = readdir(dir)) != NULL)
804         {
805             i=rindex(ent->d_name,'.');
806             if(i==NULL) continue;
807             if(strcmp(i+1,"vtu")==0){
808                     j=rindex(ent->d_name,'_');
809                     if(j==NULL) continue;
810                     number=strndup(j+1,j-i); 
a6b1b5 811                 fprintf(fh,"<DataSet timestep=\"%u\" group=\"\" part=\"0\" file=\"%s\"/>\n",atoi(number),ent->d_name);
d7639a 812                 tstep++;
SP 813                     free(number);
814             }  
815         }
816     }
f74313 817     free(dir);
d7639a 818     fprintf(fh,"</Collection>\n</VTKFile>\n");
SP 819     fclose(fh);
820     return TS_SUCCESS;
821 }
822
0a2c81 823 ts_bool write_vertex_xml_file(ts_vesicle *vesicle, ts_uint timestepno, ts_cluster_list *cstlist){
a6b1b5 824     ts_vertex_list *vlist=vesicle->vlist;
SP 825     ts_bond_list *blist=vesicle->blist;
826     ts_vertex **vtx=vlist->vtx;
40aa5b 827     ts_uint i,j;
267db5 828         char filename[10000];
SP 829         char just_name[255];
d7639a 830     FILE *fh;
267db5 831         strcpy(filename,command_line_args.path);
SP 832         sprintf(just_name,"timestep_%.6u.vtu",timestepno);
833         strcat(filename,just_name);
d7639a 834
SP 835     fh=fopen(filename, "w");
836     if(fh==NULL){
837         err("Cannot open file %s for writing");
838         return TS_FAIL;
839     }
840     /* Here comes header of the file */
40aa5b 841
SP 842     //find number of extra vtxs and bonds of polymeres
58230a 843     ts_uint monono=0, polyno=0, poly_idx=0, filno=0, fonono=0;
M 844     ts_bool poly=0, fil=0;
40aa5b 845     if(vesicle->poly_list!=NULL){
SP 846         if(vesicle->poly_list->poly[0]!=NULL){
847         polyno=vesicle->poly_list->n;
848         monono=vesicle->poly_list->poly[0]->vlist->n;
849         poly=1;
850         }
851     }
58230a 852
M 853     if(vesicle->filament_list!=NULL){
854         if(vesicle->filament_list->poly[0]!=NULL){
855         filno=vesicle->filament_list->n;
856         fonono=vesicle->filament_list->poly[0]->vlist->n;
857         fil=1;
858         }
859     }
860
854cb6 861     fprintf(fh, "<?xml version=\"1.0\"?>\n<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n");
SP 862     xml_trisurf_data(fh,vesicle);
863     fprintf(fh, " <UnstructuredGrid>\n");
8db203 864     fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+monono*polyno+fonono*filno, blist->n+monono*polyno+filno*(fonono-1)+vesicle->tlist->n);
0634b9 865 <<<<<<< HEAD
d7639a 866     fprintf(fh,"<PointData Scalars=\"scalars\">\n<DataArray type=\"Int64\" Name=\"scalars\" format=\"ascii\">");
0634b9 867 =======
5c64e2 868     fprintf(fh,"<PointData Scalars=\"vertices_idx\">\n<DataArray type=\"Int64\" Name=\"vertices_idx\" format=\"ascii\">");
0634b9 869 >>>>>>> nirgov
d7639a 870        for(i=0;i<vlist->n;i++){
a6b1b5 871         fprintf(fh,"%u ",vtx[i]->idx);
d7639a 872     }
40aa5b 873     //polymeres
SP 874     if(poly){
3c1ac1 875         poly_idx=vlist->n;
40aa5b 876         for(i=0;i<vesicle->poly_list->n;i++){
3c1ac1 877             for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){
58230a 878                 fprintf(fh,"%u ", poly_idx);
M 879             }
880         }
881     }
882     //filaments
883     if(fil){
884         poly_idx=vlist->n+monono*polyno;
885         for(i=0;i<vesicle->filament_list->n;i++){
886             for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){
887     //    fprintf(stderr,"was here\n");
3c1ac1 888                 fprintf(fh,"%u ", poly_idx);
40aa5b 889             }
SP 890         }
891     }
d7639a 892
45c708 893         fprintf(fh,"</DataArray>\n");
0634b9 894 <<<<<<< HEAD
SP 895 =======
0a2c81 896     if(cstlist!=NULL){
SP 897         fprintf(fh,"<DataArray type=\"Int64\" Name=\"vertices_in_cluster\" format=\"ascii\">");
898         for(i=0;i<vlist->n;i++){
899             if(vtx[i]->cluster!=NULL){
900                 fprintf(fh,"%u ",vtx[i]->cluster->nvtx);
901             } else {
902                 fprintf(fh,"-1 ");
903             }
904             }
905         //polymeres
906         if(poly){
907             poly_idx=vlist->n;
908             for(i=0;i<vesicle->poly_list->n;i++){
909                 for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){
910                     fprintf(fh,"-1 ");
911                 }
912             }
913         }
914         //filaments
915         if(fil){
916             poly_idx=vlist->n+monono*polyno;
917             for(i=0;i<vesicle->filament_list->n;i++){
918                 for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){
919         //    fprintf(stderr,"was here\n");
920                     fprintf(fh,"-1 ");
921                 }
922             }
923         }
924
925         fprintf(fh,"</DataArray>\n");
926
927
928     }
929
45c708 930     //here comes additional data as needed. Currently only spontaneous curvature
SP 931     fprintf(fh,"<DataArray type=\"Float64\" Name=\"spontaneous_curvature\" format=\"ascii\">");
932     for(i=0;i<vlist->n;i++){
933         fprintf(fh,"%.17e ",vtx[i]->c);
934     }
935         //polymeres
936         if(poly){
937             poly_idx=vlist->n;
938             for(i=0;i<vesicle->poly_list->n;i++){
939                 for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){
940                     fprintf(fh,"%.17e ", vesicle->poly_list->poly[i]->vlist->vtx[j]->c);
941                 }
942             }
943         }
944         //filaments
945         if(fil){
946             poly_idx=vlist->n+monono*polyno;
947             for(i=0;i<vesicle->filament_list->n;i++){
948                 for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){
949         //    fprintf(stderr,"was here\n");
950                     fprintf(fh,"%.17e ",  vesicle->filament_list->poly[i]->vlist->vtx[j]->c);
951                 }
952             }
953         }
954     fprintf(fh,"</DataArray>\n");
955
5c64e2 956     //here comes additional data. Energy!
SP 957     fprintf(fh,"<DataArray type=\"Float64\" Name=\"bending_energy\" format=\"ascii\">");
958     for(i=0;i<vlist->n;i++){
959         fprintf(fh,"%.17e ",vtx[i]->energy*vtx[i]->xk);
960     }
961         //polymeres
962         if(poly){
963             poly_idx=vlist->n;
964             for(i=0;i<vesicle->poly_list->n;i++){
965                 for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){
966                     fprintf(fh,"%.17e ", vesicle->poly_list->poly[i]->vlist->vtx[j]->energy* vesicle->poly_list->poly[i]->k);
967                 }
968             }
969         }
970         //filaments
971         if(fil){
972             poly_idx=vlist->n+monono*polyno;
973             for(i=0;i<vesicle->filament_list->n;i++){
974                 for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){
975         //    fprintf(stderr,"was here\n");
976                     fprintf(fh,"%.17e ",  vesicle->filament_list->poly[i]->vlist->vtx[j]->energy*  vesicle->filament_list->poly[i]->k);
977                 }
978             }
979         }
980     fprintf(fh,"</DataArray>\n");
981
982
0634b9 983 >>>>>>> nirgov
45c708 984     
SP 985     fprintf(fh,"</PointData>\n<CellData>\n</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n");
d7639a 986     for(i=0;i<vlist->n;i++){
aad500 987         fprintf(fh,"%.17e %.17e %.17e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
40aa5b 988     }
SP 989     //polymeres
990     if(poly){
991         for(i=0;i<vesicle->poly_list->n;i++){
992             for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
aad500 993                 fprintf(fh,"%.17e %.17e %.17e\n", vesicle->poly_list->poly[i]->vlist->vtx[j]->x,vesicle->poly_list->poly[i]->vlist->vtx[j]->y, vesicle->poly_list->poly[i]->vlist->vtx[j]->z );
58230a 994             }
M 995         }
996     }
997     //filaments
998     if(fil){
999         for(i=0;i<vesicle->filament_list->n;i++){
1000             for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){
aad500 1001                 fprintf(fh,"%.17e %.17e %.17e\n", vesicle->filament_list->poly[i]->vlist->vtx[j]->x,vesicle->filament_list->poly[i]->vlist->vtx[j]->y, vesicle->filament_list->poly[i]->vlist->vtx[j]->z );
40aa5b 1002             }
SP 1003         }
d7639a 1004     }
SP 1005
1006     fprintf(fh,"</DataArray>\n</Points>\n<Cells>\n<DataArray type=\"Int64\" Name=\"connectivity\" format=\"ascii\">");
1007     for(i=0;i<blist->n;i++){
e016c4 1008             fprintf(fh,"%u %u\n",blist->bond[i]->vtx1->idx,blist->bond[i]->vtx2->idx);
d7639a 1009     }
40aa5b 1010     //polymeres
SP 1011     if(poly){
3c1ac1 1012         poly_idx=vlist->n;
40aa5b 1013         for(i=0;i<vesicle->poly_list->n;i++){
SP 1014             for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
3c1ac1 1015 //                fprintf(fh,"%u %u\n", vesicle->poly_list->poly[i]->blist->bond[j]->vtx1->idx,vesicle->poly_list->poly[i]->blist->bond[j]->vtx2->idx);
SP 1016                 fprintf(fh,"%u %u\n", vesicle->poly_list->poly[i]->blist->bond[j]->vtx1->idx+vlist->n+i*monono,vesicle->poly_list->poly[i]->blist->bond[j]->vtx2->idx+vlist->n+i*monono);
40aa5b 1017             }
SP 1018     //grafted bonds
3c1ac1 1019         fprintf(fh,"%u %u\n", vesicle->poly_list->poly[i]->grafted_vtx->idx, vesicle->poly_list->poly[i]->vlist->vtx[0]->idx+vlist->n+i*monono);
40aa5b 1020         }
SP 1021
1022     }
1023     
58230a 1024     //filaments
M 1025     if(fil){
1026         poly_idx=vlist->n+monono*polyno;
1027         for(i=0;i<vesicle->filament_list->n;i++){
1028             for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){
1029                 fprintf(fh,"%u %u\n", vesicle->filament_list->poly[i]->blist->bond[j]->vtx1->idx+vlist->n+monono*polyno+i*fonono,vesicle->filament_list->poly[i]->blist->bond[j]->vtx2->idx+vlist->n+monono*polyno+i*fonono);
1030 //        fprintf(stderr,"was here\n");
1031             
1032             }
1033         }
1034
1035     }
8db203 1036     for(i=0;i<vesicle->tlist->n;i++){
SP 1037         fprintf(fh,"%u %u %u\n", vesicle->tlist->tria[i]->vertex[0]->idx, vesicle->tlist->tria[i]->vertex[1]->idx, vesicle->tlist->tria[i]->vertex[2]->idx);
1038     }
d7639a 1039     fprintf(fh,"</DataArray>\n<DataArray type=\"Int64\" Name=\"offsets\" format=\"ascii\">");
58230a 1040     for (i=2;i<(blist->n+monono*polyno+(fonono-1)*filno)*2+1;i+=2){
d7639a 1041     fprintf(fh,"%u ",i);
SP 1042     }
8db203 1043     for(j=i+1;j<i+3*(vesicle->tlist->n);j+=3){ //let's continue counting from where we left of
SP 1044         fprintf(fh,"%u ", j);
1045     }
d7639a 1046     fprintf(fh,"\n");
SP 1047     fprintf(fh,"</DataArray>\n<DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">\n");
8db203 1048      for (i=0;i<blist->n+monono*polyno+(fonono-1)*filno;i++){
d7639a 1049         fprintf(fh,"3 ");
SP 1050     }
8db203 1051     for(i=0;i<vesicle->tlist->n;i++){
SP 1052         fprintf(fh,"5 ");
1053     }
d7639a 1054     fprintf(fh,"</DataArray>\n</Cells>\n</Piece>\n</UnstructuredGrid>\n</VTKFile>\n");
SP 1055     fclose(fh);
1056     return TS_SUCCESS;
1057
1058 }
1059
1060 ts_bool write_vertex_vtk_file(ts_vesicle *vesicle,ts_char *filename, ts_char *text){
a6b1b5 1061     ts_vertex_list *vlist=vesicle->vlist;
SP 1062     ts_bond_list *blist=vesicle->blist;
1063     ts_vertex **vtx=vlist->vtx;
1064     ts_uint i;
d7639a 1065     FILE *fh;
SP 1066     
1067     fh=fopen(filename, "w");
1068     if(fh==NULL){
1069         err("Cannot open file %s for writing");
1070         return TS_FAIL;
1071     }
1072     /* Here comes header of the file */
1073 //    fprintf(stderr,"NSHELL=%u\n",nshell);
1074     fprintf(fh, "# vtk DataFile Version 2.0\n");
1075     /* TODO: Do a sanity check on text. Max 255 char, must not me \n terminated */ 
1076     fprintf(fh, "%s\n", text);
1077     fprintf(fh,"ASCII\n");
1078     fprintf(fh,"DATASET UNSTRUCTURED_GRID\n");
1079     fprintf(fh,"POINTS %u double\n", vlist->n);
1080     for(i=0;i<vlist->n;i++){
8f6a69 1081         fprintf(fh,"%e %e %e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
d7639a 1082     }
SP 1083     
1084     fprintf(fh,"CELLS %u %u\n",blist->n,3*blist->n);
1085     for(i=0;i<blist->n;i++){
e016c4 1086             fprintf(fh,"2 %u %u\n",blist->bond[i]->vtx1->idx,blist->bond[i]->vtx2->idx);
d7639a 1087     }
SP 1088     fprintf(fh,"CELL_TYPES %u\n",blist->n);
1089     for(i=0;i<blist->n;i++)
1090         fprintf(fh,"3\n");
1091
1092     fprintf(fh,"POINT_DATA %u\n", vlist->n);
1093     fprintf(fh,"SCALARS scalars long 1\n");
1094     fprintf(fh,"LOOKUP_TABLE default\n");
1095
1096     for(i=0;i<vlist->n;i++)
8f6a69 1097         fprintf(fh,"%u\n",vtx[i]->idx);
d7639a 1098
SP 1099     fclose(fh);
1100     return TS_SUCCESS;
1101 }
1102
1103
1104
144784 1105 ts_bool write_pov_file(ts_vesicle *vesicle, char *filename){
SP 1106     FILE *fh;
1107     ts_uint i;
1108     
1109     fh=fopen(filename, "w");
1110     if(fh==NULL){
1111         err("Cannot open file %s for writing");
1112         return TS_FAIL;
1113     }
1114
1115     for(i=0;i<vesicle->tlist->n;i++){
1116     
1117     fprintf(fh,"\ttriangle {");
1118     fprintf(fh,"\t<%e,%e,%e> <%e,%e,%e> <%e,%e,%e> }\n", 
1119     vesicle->tlist->tria[i]->vertex[0]->x,
1120     vesicle->tlist->tria[i]->vertex[0]->y,
1121     vesicle->tlist->tria[i]->vertex[0]->z,
1122
1123     vesicle->tlist->tria[i]->vertex[1]->x,
1124     vesicle->tlist->tria[i]->vertex[1]->y,
1125     vesicle->tlist->tria[i]->vertex[1]->z,
1126
1127     vesicle->tlist->tria[i]->vertex[2]->x,
1128     vesicle->tlist->tria[i]->vertex[2]->y,
1129     vesicle->tlist->tria[i]->vertex[2]->z
1130     );
1131     }
1132         
1133     fclose(fh);
1134     return TS_SUCCESS;
1135 }
1136
1137
1ab449 1138 ts_tape *parsetape(char *filename){
698ae1 1139     FILE *fd = fopen (filename, "r");
SP 1140     long length;
1141     size_t size;
1142     fseek (fd, 0, SEEK_END);
1143       length = ftell (fd);
1144     fseek (fd, 0, SEEK_SET);
1145     size=fread (tapetxt, 1, length, fd);
1146     fclose(fd);
1147     if(size);
1148     ts_tape *tape=parsetapebuffer(tapetxt);
1149     return tape;
1150 }
1151
1152 ts_tape *parsetapebuffer(char *buffer){
1ab449 1153     ts_tape *tape=(ts_tape *)calloc(1,sizeof(ts_tape));
SP 1154     tape->multiprocessing=calloc(255,sizeof(char));
40aa5b 1155
d7639a 1156     cfg_opt_t opts[] = {
1ab449 1157         CFG_SIMPLE_INT("nshell", &tape->nshell),
SP 1158         CFG_SIMPLE_INT("npoly", &tape->npoly),
1159         CFG_SIMPLE_INT("nmono", &tape->nmono),
58230a 1160     CFG_SIMPLE_INT("nfil",&tape->nfil),
M 1161     CFG_SIMPLE_INT("nfono",&tape->nfono),
e98482 1162     CFG_SIMPLE_INT("internal_poly",&tape->internal_poly),
58230a 1163     CFG_SIMPLE_INT("R_nucleus",&tape->R_nucleus),
37791b 1164     CFG_SIMPLE_FLOAT("R_nucleusX",&tape->R_nucleusX),
SP 1165     CFG_SIMPLE_FLOAT("R_nucleusY",&tape->R_nucleusY),
1166     CFG_SIMPLE_FLOAT("R_nucleusZ",&tape->R_nucleusZ),
58230a 1167     CFG_SIMPLE_FLOAT("dmax", &tape->dmax),
ea1cce 1168     CFG_SIMPLE_FLOAT("dmin_interspecies", &tape->dmin_interspecies),
1ab449 1169         CFG_SIMPLE_FLOAT("xk0",&tape->xk0),
SP 1170     CFG_SIMPLE_INT("pswitch",&tape->pswitch),
9166cb 1171     CFG_SIMPLE_INT("constvolswitch",&tape->constvolswitch),
c0ae90 1172     CFG_SIMPLE_INT("constareaswitch",&tape->constareaswitch),
43c042 1173     CFG_SIMPLE_FLOAT("constvolprecision",&tape->constvolprecision),
1ab449 1174     CFG_SIMPLE_FLOAT("pressure",&tape->pressure),
SP 1175     CFG_SIMPLE_FLOAT("k_spring",&tape->kspring),
b30f45 1176     CFG_SIMPLE_FLOAT("xi",&tape->xi),
1ab449 1177         CFG_SIMPLE_FLOAT("stepsize",&tape->stepsize),
SP 1178         CFG_SIMPLE_INT("nxmax", &tape->ncxmax),
1179         CFG_SIMPLE_INT("nymax", &tape->ncymax),
1180         CFG_SIMPLE_INT("nzmax", &tape->nczmax),
1181         CFG_SIMPLE_INT("iterations",&tape->iterations),
1182     CFG_SIMPLE_INT("mcsweeps",&tape->mcsweeps),
1183     CFG_SIMPLE_INT("inititer", &tape->inititer),
1184         CFG_SIMPLE_BOOL("quiet",&tape->quiet),
1185         CFG_SIMPLE_STR("multiprocessing",tape->multiprocessing),
1186         CFG_SIMPLE_INT("smp_cores",&tape->brezveze0),
1187         CFG_SIMPLE_INT("cluster_nodes",&tape->brezveze1),
1188         CFG_SIMPLE_INT("distributed_processes",&tape->brezveze2),
dc77e8 1189     CFG_SIMPLE_INT("spherical_harmonics_coefficients",&tape->shc),
e5858f 1190     CFG_SIMPLE_INT("number_of_vertices_with_c0", &tape->number_of_vertices_with_c0),
SP 1191     CFG_SIMPLE_FLOAT("c0",&tape->c0),
1192     CFG_SIMPLE_FLOAT("w",&tape->w),
250de4 1193     CFG_SIMPLE_FLOAT("F",&tape->F),
d7639a 1194         CFG_END()
SP 1195     };
1196     cfg_t *cfg;    
1197     ts_uint retval;
1198     cfg = cfg_init(opts, 0);
698ae1 1199     retval=cfg_parse_buf(cfg, buffer);
d7639a 1200     if(retval==CFG_FILE_ERROR){
a6b1b5 1201     fatal("No tape file.",100);
d7639a 1202     }
SP 1203     else if(retval==CFG_PARSE_ERROR){
1204     fatal("Invalid tape!",100);
1205     }
a2db52 1206
7b0c07 1207     /* here we override all values read from tape with values from commandline*/
SP 1208     getcmdline_tape(cfg,command_line_args.tape_opts);
d7639a 1209     cfg_free(cfg);
40aa5b 1210
SP 1211
1ab449 1212     /* global variables are set automatically */
SP 1213     quiet=tape->quiet;
1214     return tape;
1215 }
d7639a 1216
1ab449 1217 ts_bool tape_free(ts_tape *tape){
SP 1218     free(tape->multiprocessing);
1219     free(tape);
1220     return TS_SUCCESS;
d7639a 1221 }
f74313 1222
SP 1223
7b0c07 1224
SP 1225 ts_bool getcmdline_tape(cfg_t *cfg, char *opts){
1226
1227     char *commands, *backup, *saveptr, *saveopptr, *command, *operator[2];
1228     ts_uint i,j;
1229     commands=(char *)malloc(10000*sizeof(char));
1230     backup=commands; //since the pointer to commands will be lost, we acquire a pointer that will serve as backup.
1231     strcpy(commands,opts);
1232     for(i=0; ;i++, commands=NULL){
1233         //breaks comma separated list of commands into specific commands.
1234         command=strtok_r(commands,",",&saveptr);    
1235         if(command==NULL) break;
1236 //        fprintf(stdout,"Command %d: %s\n",i,command);    
1237         //extracts name of command and value of command into operator[2] array.
1238         for(j=0; j<2;j++,command=NULL){
1239             operator[j]=strtok_r(command,"=",&saveopptr);
1240             if(operator[j]==NULL) break;
1241 //            fprintf(stdout," ---> Operator %d: %s\n",j,operator[j]);        
1242         }
1243         //1. check: must have 2 operators.
1244         if(j!=2) fatal("Error. Command line tape options are not formatted properly",1);
1245
1246     //    cfg_setstr(cfg,operator[0],operator[1]);
1247         cmdline_to_tape(cfg,operator[0],operator[1]);
1248         //2. check: must be named properly.
1249         //3. check: must be of right format (integer, double, string, ...)
1250         
1251     }
1252     free(backup);
1253     return TS_SUCCESS;
1254 }
1255
1256
1257 ts_bool cmdline_to_tape(cfg_t *cfg, char *key, char *val){
1258
1259     cfg_opt_t *cfg_opt=cfg_getopt(cfg,key);
1260     if(cfg_opt==NULL) fatal("Commandline tape option not recognised",1); //return TS_FAIL; 
1261     switch (cfg_opt->type){
1262         case CFGT_INT:
1263             cfg_setint(cfg,key,atol(val));
1264             break;
1265         case CFGT_FLOAT:
1266             cfg_setfloat(cfg,key,atof(val));
1267             break;
1268 /*        case CFGT_BOOL:
1269             cfg_setbool(cfg,operator[0],operator[1]);
1270             break; */
1271         case CFGT_STR:
1272             cfg_setstr(cfg,key,val);
1273             break;
1274         default:
1275             break;
1276
1277     }
1278     return TS_SUCCESS;
1279 }
1280
1281
1282
f74313 1283 ts_bool read_geometry_file(char *fname, ts_vesicle *vesicle){
SP 1284     FILE *fh;
1285     ts_uint i, nvtx,nedges,ntria;
1286     ts_uint vtxi1,vtxi2;
1287     float x,y,z;
1288     ts_vertex_list *vlist;
1289     fh=fopen(fname, "r");
1290         if(fh==NULL){
1291                 err("Cannot open file for reading... Nonexistant file?");
1292                 return TS_FAIL;
1293         }
1294     ts_uint retval;
1295     retval=fscanf(fh,"%u %u %u",&nvtx, &nedges, &ntria);
1296     vesicle->vlist=init_vertex_list(nvtx);
1297     vlist=vesicle->vlist;
1298     for(i=0;i<nvtx;i++){
8f6a69 1299    //     fscanf(fh,"%F %F %F",&vlist->vtx[i]->x,&vlist->vtx[i]->y,&vlist->vtx[i]->z);
f74313 1300        retval=fscanf(fh,"%F %F %F",&x,&y,&z);
8f6a69 1301         vlist->vtx[i]->x=x;
SP 1302         vlist->vtx[i]->y=y;
1303         vlist->vtx[i]->z=z;
f74313 1304     }
SP 1305     for(i=0;i<nedges;i++){
1306         retval=fscanf(fh,"%u %u",&vtxi1,&vtxi2);
1307         bond_add(vesicle->blist,vesicle->vlist->vtx[vtxi1-1],vesicle->vlist->vtx[vtxi2-1]);
1308     }
1309     //TODO: neighbours from bonds,
1310     //TODO: triangles from neigbours
1311
1312 //    Don't need to read triangles. Already have enough data
1313     /*
1314     for(i=0;i<ntria;i++){
1315         retval=fscanf(fh,"%u %u %u", &bi1, &bi2, &bi3);
8f6a69 1316         vtxi1=vesicle->blist->vertex1->idx;
SP 1317         vtxi2=vesicle->blist->vertex1->idx;
f74313 1318         
SP 1319     }
1320     */
41a035 1321     if(retval);
f74313 1322     fclose(fh);    
SP 1323
1324
1325
1326     return TS_SUCCESS;
1327 }