trisurf-ng.git - Gitblit

Fixes error occuring when statistics.csv is empty.

Samo Penic

2016-05-17 d665c0b94c289f6c34597af3010d6b9b38a08ac4

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_nosizeof(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->cellnosizeof(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. */
	534
	535	ts_fprintf(stderr,"\n\nhere comes the help.\n\n");
	536	fatal("Ooops, read help first",1);
	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
	604
d7639a	605	ts_bool print_vertex_list(ts_vertex_list *vlist){
SP	606	ts_uint i;
	607	printf("Number of vertices: %u\n",vlist->n);
	608	for(i=0;i<vlist->n;i++){
a6b1b5	609	printf("%u: %f %f %f\n",
8f6a69	610	vlist->vtx[i]->idx,vlist->vtx[i]->x, vlist->vtx[i]->y, vlist->vtx[i]->z);
d7639a	611	}
SP	612	return TS_SUCCESS;
	613	}
	614
	615	ts_bool print_vertex_neighbours(ts_vertex_list *vlist){
	616	ts_uint i,j;
a6b1b5	617	ts_vertex **vtx=vlist->vtx;
d7639a	618	printf("Vertex id(neigh no): (neighvertex coord) (neighvertex coord) ...\n");
SP	619	for(i=0;i<vlist->n;i++){
8f6a69	620	printf("%u(%u): ",vtx[i]->idx,vtx[i]->neigh_no);
SP	621	for(j=0;j<vtx[i]->neigh_no;j++){
	622	printf("(%f,%f,%f)",vtx[i]->neigh[j]->x,
	623	vtx[i]->neigh[j]->y,vtx[i]->neigh[j]->z);
d7639a	624	}
SP	625	printf("\n");
	626	}
	627
	628	return TS_SUCCESS;
	629	}
	630
	631	ts_bool write_vertex_fcompat_file(ts_vertex_list vlist,ts_char filename){
a6b1b5	632	ts_vertex **vtx=vlist->vtx;
d7639a	633	ts_uint i;
SP	634	FILE *fh;
	635
	636	fh=fopen(filename, "w");
	637	if(fh==NULL){
	638	err("Cannot open file %s for writing");
	639	return TS_FAIL;
	640	}
	641	for(i=0;i<vlist->n;i++)
8f6a69	642	fprintf(fh," %E\t%E\t%E\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
d7639a	643
SP	644	fclose(fh);
	645	return TS_SUCCESS;
	646	}
	647
	648
	649	ts_bool fprint_vertex_list(FILE fh,ts_vertex_list vlist){
	650	ts_uint i,j;
	651	for(i=0;i<vlist->n;i++){
8f6a69	652	fprintf(fh," %.17E\t%.17E\t%.17E\t%u\n",vlist->vtx[i]->x,
SP	653	vlist->vtx[i]->y, vlist->vtx[i]->z,
	654	vlist->vtx[i]->neigh_no);
	655	for(j=0;j<vlist->vtx[i]->neigh_no;j++){
	656	fprintf(fh,"\t%u",(ts_uint)(vlist->vtx[i]->neigh[j]->idx));
a6b1b5	657	//-vlist->vtx+1));
d7639a	658	}
SP	659	fprintf(fh,"\n");
	660	}
	661	return TS_SUCCESS;
	662	}
	663
	664	ts_bool fprint_tristar(FILE fh, ts_vesicle vesicle){
	665	ts_uint i,j;
a6b1b5	666	for(i=0;i<vesicle->vlist->n;i++){
8f6a69	667	fprintf(fh,"\t%u",vesicle->vlist->vtx[i]->tristar_no);
SP	668	for(j=0;j<vesicle->vlist->vtx[i]->tristar_no;j++){
	669	fprintf(fh,"\t%u",(ts_uint)(vesicle->vlist->vtx[i]->tristar[j]->idx));//-vesicle->tlist->tria+1));
d7639a	670	}
SP	671	fprintf(fh,"\n");
	672	}
	673	return TS_SUCCESS;
	674	}
	675
	676	ts_bool fprint_triangle_list(FILE fh, ts_vesicle vesicle){
a6b1b5	677	ts_triangle_list *tlist=vesicle->tlist;
d7639a	678	ts_uint i,j;
SP	679	for(i=0;i<tlist->n;i++){
41a035	680	fprintf(fh,"\t%u",tlist->tria[i]->neigh_no);
SP	681	for(j=0;j<tlist->tria[i]->neigh_no;j++){
	682	fprintf(fh,"\t%u",(ts_uint)(tlist->tria[i]->neigh[j]->idx));//-tlist->tria+1));
d7639a	683	}
SP	684	fprintf(fh,"\n");
	685	for(j=0;j<3;j++){
41a035	686	fprintf(fh,"\t%u",(ts_uint)(tlist->tria[i]->vertex[j]->idx));//-vesicle->vlist->vtx+1));
d7639a	687	}
SP	688	fprintf(fh,"\n");
41a035	689	fprintf(fh,"%.17E\t%.17E\t%.17E\n",tlist->tria[i]->xnorm,
SP	690	tlist->tria[i]->ynorm,tlist->tria[i]->znorm);
d7639a	691	fprintf(fh,"0.00000000000000000\n0.00000000000000000\n");
SP	692	}
	693	return TS_SUCCESS;
	694	}
	695
	696	ts_bool fprint_vertex_data(FILE fh,ts_vertex_list vlist){
	697	ts_uint i,j;
	698	for(i=0;i<vlist->n;i++){
	699	fprintf(fh," %.17E\t%.17E\t%.17E\t%.17E\t%.17E\t%u\n",
8f6a69	700	vlist->vtx[i]->xk,vlist->vtx[i]->c,vlist->vtx[i]->energy,
SP	701	vlist->vtx[i]->energy_h, vlist->vtx[i]->curvature, 0);
cef6ca	702	for(j=0;j<vlist->vtx[i]->bond_no;j++){
8f6a69	703	fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length_dual);
d7639a	704	}
SP	705	fprintf(fh,"\n");
cef6ca	706	for(j=0;j<vlist->vtx[i]->bond_no;j++){
8f6a69	707	fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length);
d7639a	708	}
SP	709	fprintf(fh,"\n");
	710	}
	711	return TS_SUCCESS;
	712	}
	713
	714	ts_bool fprint_bonds(FILE fh,ts_vesicle vesicle){
	715	ts_uint i;
a6b1b5	716	for(i=0;i<vesicle->blist->n;i++){
e016c4	717	fprintf(fh,"\t%u\t%u\n",(ts_uint)(vesicle->blist->bond[i]->vtx1->idx),
a6b1b5	718	//-vesicle->vlist->vtx+1),
e016c4	719	(ts_uint)(vesicle->blist->bond[i]->vtx2->idx));
a6b1b5	720	//-vesicle->vlist.vtx+1));
d7639a	721	}
SP	722	return TS_SUCCESS;
	723	}
	724
	725
	726	ts_bool write_dout_fcompat_file(ts_vesicle vesicle, ts_char filename){
	727	FILE *fh;
	728	fh=fopen(filename, "w");
	729	if(fh==NULL){
	730	err("Cannot open file %s for writing");
	731	return TS_FAIL;
	732	}
	733	fprintf(fh,"%.17E\n%.17E\n",vesicle->stepsize,vesicle->dmax);
a6b1b5	734	fprint_vertex_list(fh,vesicle->vlist);
d7639a	735	fprint_tristar(fh,vesicle);
SP	736	fprint_triangle_list(fh,vesicle);
a6b1b5	737	fprint_vertex_data(fh,vesicle->vlist);
d7639a	738	fprint_bonds(fh,vesicle);
SP	739	fclose(fh);
	740	return TS_SUCCESS;
	741	}
	742
	743	ts_bool read_tape_fcompat_file(ts_vesicle vesicle, ts_char filename){
	744	FILE *fh;
	745	char line[255];
	746	fh=fopen(filename, "r");
	747	if(fh==NULL){
	748	err("Cannot open file for reading... Nonexistant file?");
	749	return TS_FAIL;
	750	}
a6b1b5	751	ts_uint retval=1;
d7639a	752	while(retval!=EOF){
SP	753	retval=fscanf(fh,"%s",line);
	754
	755	fprintf(stderr,"%s",line);
	756	}
	757	fclose(fh);
	758	return TS_SUCCESS;
	759	}
	760
	761	ts_bool write_master_xml_file(ts_char *filename){
	762	FILE *fh;
	763	ts_char i,j;
	764	ts_uint tstep;
	765	ts_char *number;
	766	fh=fopen(filename, "w");
	767	if(fh==NULL){
	768	err("Cannot open file %s for writing");
	769	return TS_FAIL;
	770	}
	771
	772	fprintf(fh,"<?xml version=\"1.0\"?>\n<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n<Collection>");
267db5	773	DIR *dir = opendir(command_line_args.path);
d7639a	774	if(dir){
SP	775	struct dirent *ent;
	776	tstep=0;
	777	while((ent = readdir(dir)) != NULL)
	778	{
	779	i=rindex(ent->d_name,'.');
	780	if(i==NULL) continue;
	781	if(strcmp(i+1,"vtu")==0){
	782	j=rindex(ent->d_name,'_');
	783	if(j==NULL) continue;
	784	number=strndup(j+1,j-i);
a6b1b5	785	fprintf(fh,"<DataSet timestep=\"%u\" group=\"\" part=\"0\" file=\"%s\"/>\n",atoi(number),ent->d_name);
d7639a	786	tstep++;
SP	787	free(number);
	788	}
	789	}
	790	}
f74313	791	free(dir);
d7639a	792	fprintf(fh,"</Collection>\n</VTKFile>\n");
SP	793	fclose(fh);
	794	return TS_SUCCESS;
	795	}
	796
	797	ts_bool write_vertex_xml_file(ts_vesicle *vesicle, ts_uint timestepno){
a6b1b5	798	ts_vertex_list *vlist=vesicle->vlist;
SP	799	ts_bond_list *blist=vesicle->blist;
	800	ts_vertex **vtx=vlist->vtx;
40aa5b	801	ts_uint i,j;
267db5	802	char filename[10000];
SP	803	char just_name[255];
d7639a	804	FILE *fh;
267db5	805	strcpy(filename,command_line_args.path);
SP	806	sprintf(just_name,"timestep_%.6u.vtu",timestepno);
	807	strcat(filename,just_name);
d7639a	808
SP	809	fh=fopen(filename, "w");
	810	if(fh==NULL){
	811	err("Cannot open file %s for writing");
	812	return TS_FAIL;
	813	}
	814	/* Here comes header of the file */
40aa5b	815
SP	816	//find number of extra vtxs and bonds of polymeres
58230a	817	ts_uint monono=0, polyno=0, poly_idx=0, filno=0, fonono=0;
M	818	ts_bool poly=0, fil=0;
40aa5b	819	if(vesicle->poly_list!=NULL){
SP	820	if(vesicle->poly_list->poly[0]!=NULL){
	821	polyno=vesicle->poly_list->n;
	822	monono=vesicle->poly_list->poly[0]->vlist->n;
	823	poly=1;
	824	}
	825	}
58230a	826
M	827	if(vesicle->filament_list!=NULL){
	828	if(vesicle->filament_list->poly[0]!=NULL){
	829	filno=vesicle->filament_list->n;
	830	fonono=vesicle->filament_list->poly[0]->vlist->n;
	831	fil=1;
	832	}
	833	}
	834
854cb6	835	fprintf(fh, "<?xml version=\"1.0\"?>\n<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n");
SP	836	xml_trisurf_data(fh,vesicle);
	837	fprintf(fh, " <UnstructuredGrid>\n");
58230a	838	fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+mononopolyno+fononofilno, blist->n+mononopolyno+filno(fonono-1));
d7639a	839	fprintf(fh,"<PointData Scalars=\"scalars\">\n<DataArray type=\"Int64\" Name=\"scalars\" format=\"ascii\">");
SP	840	for(i=0;i<vlist->n;i++){
a6b1b5	841	fprintf(fh,"%u ",vtx[i]->idx);
d7639a	842	}
40aa5b	843	//polymeres
SP	844	if(poly){
3c1ac1	845	poly_idx=vlist->n;
40aa5b	846	for(i=0;i<vesicle->poly_list->n;i++){
3c1ac1	847	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){
58230a	848	fprintf(fh,"%u ", poly_idx);
M	849	}
	850	}
	851	}
	852	//filaments
	853	if(fil){
	854	poly_idx=vlist->n+monono*polyno;
	855	for(i=0;i<vesicle->filament_list->n;i++){
	856	for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){
	857	// fprintf(stderr,"was here\n");
3c1ac1	858	fprintf(fh,"%u ", poly_idx);
40aa5b	859	}
SP	860	}
	861	}
d7639a	862
SP	863	fprintf(fh,"</DataArray>\n</PointData>\n<CellData>\n</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n");
	864	for(i=0;i<vlist->n;i++){
8f6a69	865	fprintf(fh,"%e %e %e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
40aa5b	866	}
SP	867	//polymeres
	868	if(poly){
	869	for(i=0;i<vesicle->poly_list->n;i++){
	870	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
	871	fprintf(fh,"%e %e %e\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	872	}
M	873	}
	874	}
	875	//filaments
	876	if(fil){
	877	for(i=0;i<vesicle->filament_list->n;i++){
	878	for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){
	879	fprintf(fh,"%e %e %e\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	880	}
SP	881	}
d7639a	882	}
SP	883
	884	fprintf(fh,"</DataArray>\n</Points>\n<Cells>\n<DataArray type=\"Int64\" Name=\"connectivity\" format=\"ascii\">");
	885	for(i=0;i<blist->n;i++){
e016c4	886	fprintf(fh,"%u %u\n",blist->bond[i]->vtx1->idx,blist->bond[i]->vtx2->idx);
d7639a	887	}
40aa5b	888	//polymeres
SP	889	if(poly){
3c1ac1	890	poly_idx=vlist->n;
40aa5b	891	for(i=0;i<vesicle->poly_list->n;i++){
SP	892	for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
3c1ac1	893	// 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	894	fprintf(fh,"%u %u\n", vesicle->poly_list->poly[i]->blist->bond[j]->vtx1->idx+vlist->n+imonono,vesicle->poly_list->poly[i]->blist->bond[j]->vtx2->idx+vlist->n+imonono);
40aa5b	895	}
SP	896	//grafted bonds
3c1ac1	897	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	898	}
SP	899
	900	}
	901
58230a	902	//filaments
M	903	if(fil){
	904	poly_idx=vlist->n+monono*polyno;
	905	for(i=0;i<vesicle->filament_list->n;i++){
	906	for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){
	907	fprintf(fh,"%u %u\n", vesicle->filament_list->poly[i]->blist->bond[j]->vtx1->idx+vlist->n+mononopolyno+ifonono,vesicle->filament_list->poly[i]->blist->bond[j]->vtx2->idx+vlist->n+mononopolyno+ifonono);
	908	// fprintf(stderr,"was here\n");
	909
	910	}
	911	}
	912
	913	}
40aa5b	914
d7639a	915	fprintf(fh,"</DataArray>\n<DataArray type=\"Int64\" Name=\"offsets\" format=\"ascii\">");
58230a	916	for (i=2;i<(blist->n+mononopolyno+(fonono-1)filno)*2+1;i+=2){
d7639a	917	fprintf(fh,"%u ",i);
SP	918	}
	919	fprintf(fh,"\n");
	920	fprintf(fh,"</DataArray>\n<DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">\n");
58230a	921	for (i=0;i<blist->n+mononopolyno+fononofilno;i++){
d7639a	922	fprintf(fh,"3 ");
SP	923	}
	924
	925	fprintf(fh,"</DataArray>\n</Cells>\n</Piece>\n</UnstructuredGrid>\n</VTKFile>\n");
	926	fclose(fh);
	927	return TS_SUCCESS;
	928
	929	}
	930
	931	ts_bool write_vertex_vtk_file(ts_vesicle vesicle,ts_char filename, ts_char *text){
a6b1b5	932	ts_vertex_list *vlist=vesicle->vlist;
SP	933	ts_bond_list *blist=vesicle->blist;
	934	ts_vertex **vtx=vlist->vtx;
	935	ts_uint i;
d7639a	936	FILE *fh;
SP	937
	938	fh=fopen(filename, "w");
	939	if(fh==NULL){
	940	err("Cannot open file %s for writing");
	941	return TS_FAIL;
	942	}
	943	/* Here comes header of the file */
	944	// fprintf(stderr,"NSHELL=%u\n",nshell);
	945	fprintf(fh, "# vtk DataFile Version 2.0\n");
	946	/* TODO: Do a sanity check on text. Max 255 char, must not me \n terminated */
	947	fprintf(fh, "%s\n", text);
	948	fprintf(fh,"ASCII\n");
	949	fprintf(fh,"DATASET UNSTRUCTURED_GRID\n");
	950	fprintf(fh,"POINTS %u double\n", vlist->n);
	951	for(i=0;i<vlist->n;i++){
8f6a69	952	fprintf(fh,"%e %e %e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
d7639a	953	}
SP	954
	955	fprintf(fh,"CELLS %u %u\n",blist->n,3*blist->n);
	956	for(i=0;i<blist->n;i++){
e016c4	957	fprintf(fh,"2 %u %u\n",blist->bond[i]->vtx1->idx,blist->bond[i]->vtx2->idx);
d7639a	958	}
SP	959	fprintf(fh,"CELL_TYPES %u\n",blist->n);
	960	for(i=0;i<blist->n;i++)
	961	fprintf(fh,"3\n");
	962
	963	fprintf(fh,"POINT_DATA %u\n", vlist->n);
	964	fprintf(fh,"SCALARS scalars long 1\n");
	965	fprintf(fh,"LOOKUP_TABLE default\n");
	966
	967	for(i=0;i<vlist->n;i++)
8f6a69	968	fprintf(fh,"%u\n",vtx[i]->idx);
d7639a	969
SP	970	fclose(fh);
	971	return TS_SUCCESS;
	972	}
	973
	974
	975
144784	976	ts_bool write_pov_file(ts_vesicle vesicle, char filename){
SP	977	FILE *fh;
	978	ts_uint i;
	979
	980	fh=fopen(filename, "w");
	981	if(fh==NULL){
	982	err("Cannot open file %s for writing");
	983	return TS_FAIL;
	984	}
	985
	986	for(i=0;i<vesicle->tlist->n;i++){
	987
	988	fprintf(fh,"\ttriangle {");
	989	fprintf(fh,"\t<%e,%e,%e> <%e,%e,%e> <%e,%e,%e> }\n",
	990	vesicle->tlist->tria[i]->vertex[0]->x,
	991	vesicle->tlist->tria[i]->vertex[0]->y,
	992	vesicle->tlist->tria[i]->vertex[0]->z,
	993
	994	vesicle->tlist->tria[i]->vertex[1]->x,
	995	vesicle->tlist->tria[i]->vertex[1]->y,
	996	vesicle->tlist->tria[i]->vertex[1]->z,
	997
	998	vesicle->tlist->tria[i]->vertex[2]->x,
	999	vesicle->tlist->tria[i]->vertex[2]->y,
	1000	vesicle->tlist->tria[i]->vertex[2]->z
	1001	);
	1002	}
	1003
	1004	fclose(fh);
	1005	return TS_SUCCESS;
	1006	}
	1007
	1008
1ab449	1009	ts_tape parsetape(char filename){
698ae1	1010	FILE *fd = fopen (filename, "r");
SP	1011	long length;
	1012	size_t size;
	1013	fseek (fd, 0, SEEK_END);
	1014	length = ftell (fd);
	1015	fseek (fd, 0, SEEK_SET);
	1016	size=fread (tapetxt, 1, length, fd);
	1017	fclose(fd);
	1018	if(size);
	1019	ts_tape *tape=parsetapebuffer(tapetxt);
	1020	return tape;
	1021	}
	1022
	1023	ts_tape parsetapebuffer(char buffer){
1ab449	1024	ts_tape tape=(ts_tape )calloc(1,sizeof(ts_tape));
SP	1025	tape->multiprocessing=calloc(255,sizeof(char));
40aa5b	1026
d7639a	1027	cfg_opt_t opts[] = {
1ab449	1028	CFG_SIMPLE_INT("nshell", &tape->nshell),
SP	1029	CFG_SIMPLE_INT("npoly", &tape->npoly),
	1030	CFG_SIMPLE_INT("nmono", &tape->nmono),
58230a	1031	CFG_SIMPLE_INT("nfil",&tape->nfil),
M	1032	CFG_SIMPLE_INT("nfono",&tape->nfono),
	1033	CFG_SIMPLE_INT("R_nucleus",&tape->R_nucleus),
	1034	CFG_SIMPLE_FLOAT("dmax", &tape->dmax),
ea1cce	1035	CFG_SIMPLE_FLOAT("dmin_interspecies", &tape->dmin_interspecies),
1ab449	1036	CFG_SIMPLE_FLOAT("xk0",&tape->xk0),
SP	1037	CFG_SIMPLE_INT("pswitch",&tape->pswitch),
9166cb	1038	CFG_SIMPLE_INT("constvolswitch",&tape->constvolswitch),
c0ae90	1039	CFG_SIMPLE_INT("constareaswitch",&tape->constareaswitch),
43c042	1040	CFG_SIMPLE_FLOAT("constvolprecision",&tape->constvolprecision),
1ab449	1041	CFG_SIMPLE_FLOAT("pressure",&tape->pressure),
SP	1042	CFG_SIMPLE_FLOAT("k_spring",&tape->kspring),
b30f45	1043	CFG_SIMPLE_FLOAT("xi",&tape->xi),
1ab449	1044	CFG_SIMPLE_FLOAT("stepsize",&tape->stepsize),
SP	1045	CFG_SIMPLE_INT("nxmax", &tape->ncxmax),
	1046	CFG_SIMPLE_INT("nymax", &tape->ncymax),
	1047	CFG_SIMPLE_INT("nzmax", &tape->nczmax),
	1048	CFG_SIMPLE_INT("iterations",&tape->iterations),
	1049	CFG_SIMPLE_INT("mcsweeps",&tape->mcsweeps),
	1050	CFG_SIMPLE_INT("inititer", &tape->inititer),
	1051	CFG_SIMPLE_BOOL("quiet",&tape->quiet),
	1052	CFG_SIMPLE_STR("multiprocessing",tape->multiprocessing),
	1053	CFG_SIMPLE_INT("smp_cores",&tape->brezveze0),
	1054	CFG_SIMPLE_INT("cluster_nodes",&tape->brezveze1),
	1055	CFG_SIMPLE_INT("distributed_processes",&tape->brezveze2),
dc77e8	1056	CFG_SIMPLE_INT("spherical_harmonics_coefficients",&tape->shc),
d7639a	1057	CFG_END()
SP	1058	};
	1059	cfg_t *cfg;
	1060	ts_uint retval;
	1061	cfg = cfg_init(opts, 0);
698ae1	1062	retval=cfg_parse_buf(cfg, buffer);
d7639a	1063	if(retval==CFG_FILE_ERROR){
a6b1b5	1064	fatal("No tape file.",100);
d7639a	1065	}
SP	1066	else if(retval==CFG_PARSE_ERROR){
	1067	fatal("Invalid tape!",100);
	1068	}
a2db52	1069
7b0c07	1070	/* here we override all values read from tape with values from commandline*/
SP	1071	getcmdline_tape(cfg,command_line_args.tape_opts);
d7639a	1072	cfg_free(cfg);
40aa5b	1073
SP	1074
1ab449	1075	/* global variables are set automatically */
SP	1076	quiet=tape->quiet;
	1077	return tape;
	1078	}
d7639a	1079
1ab449	1080	ts_bool tape_free(ts_tape *tape){
SP	1081	free(tape->multiprocessing);
	1082	free(tape);
	1083	return TS_SUCCESS;
d7639a	1084	}
f74313	1085
SP	1086
7b0c07	1087
SP	1088	ts_bool getcmdline_tape(cfg_t cfg, char opts){
	1089
	1090	char commands, backup, saveptr, saveopptr, command, operator[2];
	1091	ts_uint i,j;
	1092	commands=(char )malloc(10000sizeof(char));
	1093	backup=commands; //since the pointer to commands will be lost, we acquire a pointer that will serve as backup.
	1094	strcpy(commands,opts);
	1095	for(i=0; ;i++, commands=NULL){
	1096	//breaks comma separated list of commands into specific commands.
	1097	command=strtok_r(commands,",",&saveptr);
	1098	if(command==NULL) break;
	1099	// fprintf(stdout,"Command %d: %s\n",i,command);
	1100	//extracts name of command and value of command into operator[2] array.
	1101	for(j=0; j<2;j++,command=NULL){
	1102	operator[j]=strtok_r(command,"=",&saveopptr);
	1103	if(operator[j]==NULL) break;
	1104	// fprintf(stdout," ---> Operator %d: %s\n",j,operator[j]);
	1105	}
	1106	//1. check: must have 2 operators.
	1107	if(j!=2) fatal("Error. Command line tape options are not formatted properly",1);
	1108
	1109	// cfg_setstr(cfg,operator[0],operator[1]);
	1110	cmdline_to_tape(cfg,operator[0],operator[1]);
	1111	//2. check: must be named properly.
	1112	//3. check: must be of right format (integer, double, string, ...)
	1113
	1114	}
	1115	free(backup);
	1116	return TS_SUCCESS;
	1117	}
	1118
	1119
	1120	ts_bool cmdline_to_tape(cfg_t cfg, char key, char *val){
	1121
	1122	cfg_opt_t *cfg_opt=cfg_getopt(cfg,key);
	1123	if(cfg_opt==NULL) fatal("Commandline tape option not recognised",1); //return TS_FAIL;
	1124	switch (cfg_opt->type){
	1125	case CFGT_INT:
	1126	cfg_setint(cfg,key,atol(val));
	1127	break;
	1128	case CFGT_FLOAT:
	1129	cfg_setfloat(cfg,key,atof(val));
	1130	break;
	1131	/* case CFGT_BOOL:
	1132	cfg_setbool(cfg,operator[0],operator[1]);
	1133	break; */
	1134	case CFGT_STR:
	1135	cfg_setstr(cfg,key,val);
	1136	break;
	1137	default:
	1138	break;
	1139
	1140	}
	1141	return TS_SUCCESS;
	1142	}
	1143
	1144
	1145
f74313	1146	ts_bool read_geometry_file(char fname, ts_vesicle vesicle){
SP	1147	FILE *fh;
	1148	ts_uint i, nvtx,nedges,ntria;
	1149	ts_uint vtxi1,vtxi2;
	1150	float x,y,z;
	1151	ts_vertex_list *vlist;
	1152	fh=fopen(fname, "r");
	1153	if(fh==NULL){
	1154	err("Cannot open file for reading... Nonexistant file?");
	1155	return TS_FAIL;
	1156	}
	1157	ts_uint retval;
	1158	retval=fscanf(fh,"%u %u %u",&nvtx, &nedges, &ntria);
	1159	vesicle->vlist=init_vertex_list(nvtx);
	1160	vlist=vesicle->vlist;
	1161	for(i=0;i<nvtx;i++){
8f6a69	1162	// fscanf(fh,"%F %F %F",&vlist->vtx[i]->x,&vlist->vtx[i]->y,&vlist->vtx[i]->z);
f74313	1163	retval=fscanf(fh,"%F %F %F",&x,&y,&z);
8f6a69	1164	vlist->vtx[i]->x=x;
SP	1165	vlist->vtx[i]->y=y;
	1166	vlist->vtx[i]->z=z;
f74313	1167	}
SP	1168	for(i=0;i<nedges;i++){
	1169	retval=fscanf(fh,"%u %u",&vtxi1,&vtxi2);
	1170	bond_add(vesicle->blist,vesicle->vlist->vtx[vtxi1-1],vesicle->vlist->vtx[vtxi2-1]);
	1171	}
	1172	//TODO: neighbours from bonds,
	1173	//TODO: triangles from neigbours
	1174
	1175	// Don't need to read triangles. Already have enough data
	1176	/*
	1177	for(i=0;i<ntria;i++){
	1178	retval=fscanf(fh,"%u %u %u", &bi1, &bi2, &bi3);
8f6a69	1179	vtxi1=vesicle->blist->vertex1->idx;
SP	1180	vtxi2=vesicle->blist->vertex1->idx;
f74313	1181
SP	1182	}
	1183	*/
41a035	1184	if(retval);
f74313	1185	fclose(fh);
SP	1186
	1187
	1188
	1189	return TS_SUCCESS;
	1190	}