trisurf-ng.git - Gitblit

Samo Penic

2014-03-06 083e03adeaf985c9da1ed31d259ced9bb54c207c

commit \| author \| age
f74313	1	#include "general.h"
d7639a	2	#include<stdio.h>
SP	3	#include "io.h"
f74313	4	#include <confuse.h>
SP	5	#include "vertex.h"
	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"
d7639a	13
083e03	14	#include <getopt.h>
SP	15	#include <sys/stat.h>
	16	#include <sys/types.h>
	17	#include <dirent.h>
	18	#include <errno.h>
e9eab4	19	/ DUMP STATE TO DISK DRIVE /
SP	20
	21	ts_bool dump_state(ts_vesicle *vesicle){
	22
	23	/* save current state with wrong pointers. Will fix that later */
	24	ts_uint i,j,k;
	25	FILE *fh=fopen("dump.bin","wb");
	26
	27	/* dump vesicle */
86f5e7	28	fwrite(vesicle, sizeof(ts_vesicle),1,fh);
e9eab4	29	/* dump vertex list */
SP	30	fwrite(vesicle->vlist, sizeof(ts_vertex_list),1,fh);
	31	/* dump bond list */
	32	fwrite(vesicle->blist, sizeof(ts_bond_list),1,fh);
	33	/* dump triangle list */
	34	fwrite(vesicle->tlist, sizeof(ts_triangle_list),1,fh);
	35	/* dump cell list */
	36	fwrite(vesicle->clist, sizeof(ts_cell_list),1,fh);
	37	/* dump poly list */
	38	fwrite(vesicle->poly_list, sizeof(ts_poly_list),1,fh);
	39	/* level 1 complete */
	40
	41	/dump vertices/
	42	for(i=0;i<vesicle->vlist->n;i++){
	43	fwrite(vesicle->vlist->vtx[i],sizeof(ts_vertex),1,fh);
	44	/* dump pointer offsets for:
	45	neigh
	46	bond
	47	tria
	48	cell is ignored
	49	*/
	50	for(j=0;j<vesicle->vlist->vtx[i]->neigh_no;j++){
3c1ac1	51	fwrite(&vesicle->vlist->vtx[i]->neigh[j]->idx,sizeof(ts_uint),1,fh);
e9eab4	52	}
SP	53	for(j=0;j<vesicle->vlist->vtx[i]->bond_no;j++){
3c1ac1	54	fwrite(&vesicle->vlist->vtx[i]->bond[j]->idx,sizeof(ts_uint),1,fh);
e9eab4	55	}
SP	56	for(j=0;j<vesicle->vlist->vtx[i]->tristar_no;j++){
3c1ac1	57	fwrite(&vesicle->vlist->vtx[i]->tristar[j]->idx,sizeof(ts_uint),1,fh);
e9eab4	58	}
SP	59	}
	60
	61	/dump bonds/
	62	for(i=0;i<vesicle->blist->n;i++){
	63	fwrite(vesicle->blist->bond[i],sizeof(ts_bond),1,fh);
	64	/* dump pointer offsets for vtx1 and vtx2 */
3c1ac1	65	//off=(ts_ulong)(vesicle->blist->bond[i]->vtx1-vesicle->vlist->vtx[0]);
SP	66	fwrite(&vesicle->blist->bond[i]->vtx1->idx,sizeof(ts_uint),1,fh);
	67	//off=(ts_ulong)(vesicle->blist->bond[i]->vtx2-vesicle->vlist->vtx[0]);
	68	fwrite(&vesicle->blist->bond[i]->vtx2->idx,sizeof(ts_uint),1,fh);
e9eab4	69	}
SP	70
	71	/dump triangles/
	72	for(i=0;i<vesicle->tlist->n;i++){
	73	fwrite(vesicle->tlist->tria[i],sizeof(ts_triangle),1,fh);
	74	/* dump pointer offsets for vertex */
3c1ac1	75	fwrite(&vesicle->tlist->tria[i]->vertex[0]->idx,sizeof(ts_uint),1,fh);
SP	76	fwrite(&vesicle->tlist->tria[i]->vertex[1]->idx,sizeof(ts_uint),1,fh);
	77	fwrite(&vesicle->tlist->tria[i]->vertex[2]->idx,sizeof(ts_uint),1,fh);
e9eab4	78	/* dump pointer offsets for neigh */
SP	79	for(j=0;j<vesicle->tlist->tria[i]->neigh_no;j++){
3c1ac1	80	fwrite(&vesicle->tlist->tria[i]->neigh[j]->idx,sizeof(ts_uint),1,fh);
e9eab4	81	}
SP	82	}
	83
	84
	85	/dump polymeres /
	86	for(i=0;i<vesicle->poly_list->n;i++){
	87	fwrite(vesicle->poly_list->poly[i],sizeof(ts_poly),1,fh);
3c1ac1	88	fwrite(vesicle->poly_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh);
SP	89	fwrite(vesicle->poly_list->poly[i]->blist,sizeof(ts_bond_list),1,fh);
e9eab4	90	}
SP	91
	92	/* dump poly vertex(monomer) list*/
	93	for(i=0;i<vesicle->poly_list->n;i++){
	94	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
	95	fwrite(vesicle->poly_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh);
	96	/* dump offset for neigh and bond */
	97	for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh_no;k++){
3c1ac1	98	// off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]-vesicle->poly_list->poly[i]->vlist->vtx[0]);
SP	99	fwrite(&vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]->idx,sizeof(ts_uint),1,fh);
e9eab4	100	}
SP	101	for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->bond_no;k++){
3c1ac1	102	//off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]-vesicle->poly_list->poly[i]->blist->bond[0]);
SP	103	fwrite(&vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]->idx,sizeof(ts_uint),1,fh);
e9eab4	104	}
SP	105	}
3c1ac1	106	// grafted vtx on vesicle data dump
SP	107	fwrite(&vesicle->poly_list->poly[i]->grafted_vtx->idx, sizeof(ts_uint),1,fh);
e9eab4	108	}
SP	109	/* dump poly bonds between monomers list*/
	110	for(i=0;i<vesicle->poly_list->n;i++){
	111	for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
	112	fwrite(vesicle->poly_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh);
	113	/* dump vtx1 and vtx2 offsets */
3c1ac1	114	//off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx1-vesicle->poly_list->poly[i]->vlist->vtx[0]);
SP	115	fwrite(&vesicle->poly_list->poly[i]->blist->bond[j]->vtx1->idx,sizeof(ts_uint),1,fh);
	116	// off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx2-vesicle->poly_list->poly[i]->vlist->vtx[0]);
	117	fwrite(&vesicle->poly_list->poly[i]->blist->bond[j]->vtx2->idx,sizeof(ts_uint),1,fh);
e9eab4	118	}
SP	119	}
	120
	121	/* pointer offsets for fixing the restored pointers */
	122	/* need pointers for
	123	vlist->vtx
	124	blist->bond
	125	tlist->tria
	126	clist->cell
	127	poly_list->poly
	128	and for each poly:
	129	poly_list->poly->vtx
	130	poly_list->poly->bond
	131	*/
	132
063289	133	fwrite(vesicle->clist, sizeof(ts_cell_list),1, fh);
SP	134
e9eab4	135	fclose(fh);
SP	136	return TS_SUCCESS;
	137	}
	138
	139
	140	/ RESTORE DUMP FROM DISK /
	141	ts_vesicle *restore_state(){
	142	ts_uint i,j,k;
	143	FILE *fh=fopen("dump.bin","rb");
	144	ts_uint retval;
3c1ac1	145	ts_uint idx;
e9eab4	146
SP	147	/* we restore all the data from the dump */
	148	/* restore vesicle */
	149	ts_vesicle vesicle=(ts_vesicle )calloc(1,sizeof(ts_vesicle));
86f5e7	150	retval=fread(vesicle, sizeof(ts_vesicle),1,fh);
SP	151	// fprintf(stderr,"was here! %e\n",vesicle->dmax);
	152
e9eab4	153	/* restore vertex list */
SP	154	vesicle->vlist=(ts_vertex_list *)malloc(sizeof(ts_vertex_list));
	155	retval=fread(vesicle->vlist, sizeof(ts_vertex_list),1,fh);
	156	/* restore bond list */
	157	vesicle->blist=(ts_bond_list *)malloc(sizeof(ts_bond_list));
	158	retval=fread(vesicle->blist, sizeof(ts_bond_list),1,fh);
	159	/* restore triangle list */
	160	vesicle->tlist=(ts_triangle_list *)malloc(sizeof(ts_triangle_list));
	161	retval=fread(vesicle->tlist, sizeof(ts_triangle_list),1,fh);
	162	/* restore cell list */
	163	vesicle->clist=(ts_cell_list *)malloc(sizeof(ts_cell_list));
	164	retval=fread(vesicle->clist, sizeof(ts_cell_list),1,fh);
	165	/* restore poly list */
	166	vesicle->poly_list=(ts_poly_list *)calloc(1,sizeof(ts_poly_list));
	167	retval=fread(vesicle->poly_list, sizeof(ts_poly_list),1,fh);
	168	/* level 1 complete */
	169
3c1ac1	170	/* prerequisity. Bonds must be malloced before vertexes are recreated */
SP	171	vesicle->blist->bond=(ts_bond *)calloc(vesicle->blist->n,sizeof(ts_bond ));
	172	for(i=0;i<vesicle->blist->n;i++){
	173	vesicle->blist->bond[i]=(ts_bond *)malloc(sizeof(ts_bond));
	174	}
	175	/* prerequisity. Triangles must be malloced before vertexes are recreated */
	176	vesicle->tlist->tria=(ts_triangle *)calloc(vesicle->tlist->n,sizeof(ts_triangle ));
	177	for(i=0;i<vesicle->tlist->n;i++){
	178	vesicle->tlist->tria[i]=(ts_triangle *)malloc(sizeof(ts_triangle));
	179	}
a37ba2	180	/* prerequisity. Vertices must be malloced before vertexes are recreated */
e9eab4	181	vesicle->vlist->vtx=(ts_vertex *)calloc(vesicle->vlist->n,sizeof(ts_vertex ));
a37ba2	182	for(i=0;i<vesicle->vlist->n;i++){
e9eab4	183	vesicle->vlist->vtx[i]=(ts_vertex *)malloc(sizeof(ts_vertex));
a37ba2	184	}
SP	185	/restore vertices/
	186	for(i=0;i<vesicle->vlist->n;i++){
e9eab4	187	retval=fread(vesicle->vlist->vtx[i],sizeof(ts_vertex),1,fh);
SP	188	/restore neigh, bond, tristar. Ignoring cell /
	189	vesicle->vlist->vtx[i]->neigh=(ts_vertex *)calloc(vesicle->vlist->vtx[i]->neigh_no, sizeof(ts_vertex ));
	190	for(j=0;j<vesicle->vlist->vtx[i]->neigh_no;j++){
3c1ac1	191	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	192	vesicle->vlist->vtx[i]->neigh[j]=vesicle->vlist->vtx[idx];
e9eab4	193	}
SP	194	vesicle->vlist->vtx[i]->bond=(ts_bond *)calloc(vesicle->vlist->vtx[i]->bond_no, sizeof(ts_bond ));
	195	for(j=0;j<vesicle->vlist->vtx[i]->bond_no;j++){
3c1ac1	196	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	197	/* 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 */
	198	vesicle->vlist->vtx[i]->bond[j]=vesicle->blist->bond[idx];
e9eab4	199	}
SP	200
	201	vesicle->vlist->vtx[i]->tristar=(ts_triangle *)calloc(vesicle->vlist->vtx[i]->tristar_no, sizeof(ts_triangle ));
	202	for(j=0;j<vesicle->vlist->vtx[i]->tristar_no;j++){
3c1ac1	203	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	204	/* same comment as above */
	205	vesicle->vlist->vtx[i]->tristar[j]=vesicle->tlist->tria[idx];
e9eab4	206	}
SP	207
	208	}
	209
	210	/restore bonds/
3c1ac1	211	// vesicle->blist->bond=(ts_bond *)calloc(vesicle->blist->n,sizeof(ts_bond )); // done before.
e9eab4	212	for(i=0;i<vesicle->blist->n;i++){
3c1ac1	213	// vesicle->blist->bond[i]=(ts_bond *)malloc(sizeof(ts_bond)); //done before.
e9eab4	214	retval=fread(vesicle->blist->bond[i],sizeof(ts_bond),1,fh);
SP	215	/* restore vtx1 and vtx2 */
3c1ac1	216	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	217	vesicle->blist->bond[i]->vtx1=vesicle->vlist->vtx[idx];
	218	retval=fread(&idx,sizeof(ts_uint),1,fh);
	219	vesicle->blist->bond[i]->vtx2=vesicle->vlist->vtx[idx];
e9eab4	220	}
SP	221
	222	/restore triangles/
3c1ac1	223	// vesicle->tlist->tria=(ts_triangle *)calloc(vesicle->tlist->n,sizeof(ts_triangle )); // done before
e9eab4	224	for(i=0;i<vesicle->tlist->n;i++){
3c1ac1	225	// vesicle->tlist->tria[i]=(ts_triangle *)malloc(sizeof(ts_triangle)); // done before
e9eab4	226	retval=fread(vesicle->tlist->tria[i],sizeof(ts_triangle),1,fh);
SP	227	/* restore pointers for vertices */
3c1ac1	228	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	229	vesicle->tlist->tria[i]->vertex[0]=vesicle->vlist->vtx[idx];
	230	retval=fread(&idx,sizeof(ts_uint),1,fh);
	231	vesicle->tlist->tria[i]->vertex[1]=vesicle->vlist->vtx[idx];
	232	retval=fread(&idx,sizeof(ts_uint),1,fh);
	233	vesicle->tlist->tria[i]->vertex[2]=vesicle->vlist->vtx[idx];
e9eab4	234	/* restore pointers for neigh */
3c1ac1	235	vesicle->tlist->tria[i]->neigh=(ts_triangle *)malloc(vesicle->tlist->tria[i]->neigh_nosizeof(ts_triangle *));
e9eab4	236	for(j=0;j<vesicle->tlist->tria[i]->neigh_no;j++){
3c1ac1	237	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	238	vesicle->tlist->tria[i]->neigh[j]=vesicle->tlist->tria[idx];
e9eab4	239	}
SP	240
	241	}
	242
	243	/restore cells /
	244	/TODO: do we need to recalculate cells here? /
	245	/* vesicle->clist->cell=(ts_cell *)malloc(vesicle->clist->cellnosizeof(ts_cell *));
	246	for(i=0;i<vesicle->clist->cellno;i++){
	247	vesicle->clist->cell[i]=(ts_cell *)malloc(sizeof(ts_cell));
	248	retval=fread(vesicle->clist->cell[i],sizeof(ts_cell),1,fh);
	249	}
	250	*/
	251	/restore polymeres /
	252	vesicle->poly_list->poly = (ts_poly *)calloc(vesicle->poly_list->n,sizeof(ts_poly ));
	253	for(i=0;i<vesicle->poly_list->n;i++){
	254	vesicle->poly_list->poly[i]=(ts_poly *)calloc(1,sizeof(ts_poly));
	255	retval=fread(vesicle->poly_list->poly[i],sizeof(ts_poly),1,fh);
3c1ac1	256	vesicle->poly_list->poly[i]->vlist=(ts_vertex_list *)calloc(1,sizeof(ts_vertex_list));
SP	257	retval=fread(vesicle->poly_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh);
	258	vesicle->poly_list->poly[i]->blist=(ts_bond_list *)calloc(1,sizeof(ts_bond_list));
	259	retval=fread(vesicle->poly_list->poly[i]->blist,sizeof(ts_bond_list),1,fh);
	260	/* initialize adress space for pointers that will hold specific vertices (monomers) and bonds */
	261	vesicle->poly_list->poly[i]->vlist->vtx=(ts_vertex *)calloc(vesicle->poly_list->poly[i]->vlist->n,sizeof(ts_vertex ));
	262	vesicle->poly_list->poly[i]->blist->bond=(ts_bond *)calloc(vesicle->poly_list->poly[i]->blist->n,sizeof(ts_bond ));
	263	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
	264	vesicle->poly_list->poly[i]->vlist->vtx[j]=(ts_vertex *)malloc(sizeof(ts_vertex));
	265	}
	266	for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
	267	vesicle->poly_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond));
	268	}
	269
e9eab4	270	}
3c1ac1	271
e9eab4	272
SP	273	/* restore poly vertex(monomer) list*/
	274	for(i=0;i<vesicle->poly_list->n;i++){
	275	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
	276	retval=fread(vesicle->poly_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh);
3c1ac1	277
e9eab4	278	/* restore neigh and bonds */
SP	279	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 ));
	280	for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh_no;k++){
3c1ac1	281	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	282	vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]=vesicle->poly_list->poly[i]->vlist->vtx[idx];
e9eab4	283	}
SP	284	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 ));
	285	for(k=0;k<vesicle->poly_list->poly[i]->vlist->vtx[j]->bond_no;k++){
3c1ac1	286	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	287	vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]=vesicle->poly_list->poly[i]->blist->bond[idx];
e9eab4	288	}
SP	289
	290	}
3c1ac1	291	/* restore grafted vtx on vesicle and grafted_poly */
SP	292	retval=fread(&idx,sizeof(ts_uint),1,fh);
	293	vesicle->vlist->vtx[idx]->grafted_poly=vesicle->poly_list->poly[i];
	294	vesicle->poly_list->poly[i]->grafted_vtx=vesicle->vlist->vtx[idx];
e9eab4	295	}
3c1ac1	296
e9eab4	297	/* restore poly bonds between monomers list*/
SP	298	for(i=0;i<vesicle->poly_list->n;i++){
	299	for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
3c1ac1	300	// vesicle->poly_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond));
e9eab4	301	retval=fread(vesicle->poly_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh);
SP	302	/* restore vtx1 and vtx2 */
3c1ac1	303	retval=fread(&idx,sizeof(ts_uint),1,fh);
SP	304	vesicle->poly_list->poly[i]->blist->bond[j]->vtx1=vesicle->poly_list->poly[i]->vlist->vtx[idx];
	305	retval=fread(&idx,sizeof(ts_uint),1,fh);
	306	vesicle->poly_list->poly[i]->blist->bond[j]->vtx2=vesicle->poly_list->poly[i]->vlist->vtx[idx];
e9eab4	307	}
SP	308	}
	309
063289	310	// recreating space for cells //
SP	311	vesicle->clist=(ts_cell_list *)malloc(sizeof(ts_cell_list));
	312	retval=fread(vesicle->clist, sizeof(ts_cell_list), 1,fh);
	313	vesicle->clist->cell=(ts_cell *)malloc(sizeof(ts_cell )vesicle->clist->ncmax[0]vesicle->clist->ncmax[1]*vesicle->clist->ncmax[2]);
	314	for(i=0;i<vesicle->clist->ncmax[0]vesicle->clist->ncmax[1]vesicle->clist->ncmax[2];i++){
	315	vesicle->clist->cell[i]=(ts_cell *)calloc(1,sizeof(ts_cell));
	316	vesicle->clist->cell[i]->idx=i+1; // We enumerate cells! Probably never required!
	317	}
e9eab4	318
SP	319	if(retval);
	320	fclose(fh);
	321	return vesicle;
	322	}
	323
	324
	325
083e03	326	ts_bool parse_args(int argc, char **argv){
SP	327	int c, retval;
	328	DIR *dir;
	329	path[0]=0;
	330	while (1)
	331	{
	332	static struct option long_options[] =
	333	{
	334	{"force-from-tape", no_argument, &force_from_tape, 1},
	335	{"tape", no_argument, 0, 't'},
	336	{"output-file", required_argument, 0, 'o'},
	337	{"directory", required_argument, 0, 'd'},
	338	{"dump-file", required_argument,0, 'f'},
	339	{0, 0, 0, 0}
	340	};
	341	/* getopt_long stores the option index here. */
	342	int option_index = 0;
	343
	344	c = getopt_long (argc, argv, "d:fot",
	345	long_options, &option_index);
	346
	347	/* Detect the end of the options. */
	348	if (c == -1)
	349	break;
	350
	351	switch (c)
	352	{
	353	case 0:
	354	/* If this option set a flag, do nothing else now. */
	355	if (long_options[option_index].flag != 0)
	356	break;
	357	printf ("option %s", long_options[option_index].name);
	358	if (optarg)
	359	printf (" with arg %s", optarg);
	360	printf ("\n");
	361	break;
	362
	363	case 't':
	364	//check if tape exists. If not, fail immediately.
	365	puts ("option -t\n");
	366	break;
	367
	368	case 'o':
	369	//set filename of master output file
	370	printf ("option -o with value `%s'\n", optarg);
	371	break;
	372
	373	case 'd':
	374	//check if directory exists. If not create one. If creation is
	375	//successful, set directory for output files.
	376	//printf ("option -d with value `%s'\n", optarg);
	377	dir = opendir(optarg);
	378	if (dir)
	379	{
	380	/* Directory exists. */
	381	closedir(dir);
	382	}
	383	else if (ENOENT == errno)
	384	{
	385	/* Directory does not exist. */
	386	retval=mkdir(optarg, 0700);
	387	if(retval){
	388	fatal("Could not create requested directory. Check if you have permissions",1);
	389	}
	390	}
	391	else
	392	{
	393	/* opendir() failed for some other reason. */
	394	fatal("Could not check if directory exists. Reason unknown",1);
	395	}
	396	ts_fprintf(stdout,"\n*** Using output directory: %s\n\n", optarg);
	397	// sprintf(path,"%s", optarg);
	398	strcpy(path, optarg);
	399	// ts_fprintf(stdout,"ok!\n");
	400
	401	break;
	402
	403	case 'f':
	404	//check if dump file specified exists. Defaults to dump.bin
	405	break;
	406
	407	case '?':
	408	/* getopt_long already printed an error message. */
	409
	410	ts_fprintf(stderr,"\n\nhere comes the help.\n\n");
	411	fatal("Ooops, read help first",1);
	412	break;
	413
	414	default:
	415	exit (1);
	416	}
	417	}
	418
	419	return TS_SUCCESS;
	420
	421	}
	422
	423
	424
	425
d7639a	426	ts_bool print_vertex_list(ts_vertex_list *vlist){
SP	427	ts_uint i;
	428	printf("Number of vertices: %u\n",vlist->n);
	429	for(i=0;i<vlist->n;i++){
a6b1b5	430	printf("%u: %f %f %f\n",
8f6a69	431	vlist->vtx[i]->idx,vlist->vtx[i]->x, vlist->vtx[i]->y, vlist->vtx[i]->z);
d7639a	432	}
SP	433	return TS_SUCCESS;
	434	}
	435
	436	ts_bool print_vertex_neighbours(ts_vertex_list *vlist){
	437	ts_uint i,j;
a6b1b5	438	ts_vertex **vtx=vlist->vtx;
d7639a	439	printf("Vertex id(neigh no): (neighvertex coord) (neighvertex coord) ...\n");
SP	440	for(i=0;i<vlist->n;i++){
8f6a69	441	printf("%u(%u): ",vtx[i]->idx,vtx[i]->neigh_no);
SP	442	for(j=0;j<vtx[i]->neigh_no;j++){
	443	printf("(%f,%f,%f)",vtx[i]->neigh[j]->x,
	444	vtx[i]->neigh[j]->y,vtx[i]->neigh[j]->z);
d7639a	445	}
SP	446	printf("\n");
	447	}
	448
	449	return TS_SUCCESS;
	450	}
	451
	452	ts_bool write_vertex_fcompat_file(ts_vertex_list vlist,ts_char filename){
a6b1b5	453	ts_vertex **vtx=vlist->vtx;
d7639a	454	ts_uint i;
SP	455	FILE *fh;
	456
	457	fh=fopen(filename, "w");
	458	if(fh==NULL){
	459	err("Cannot open file %s for writing");
	460	return TS_FAIL;
	461	}
	462	for(i=0;i<vlist->n;i++)
8f6a69	463	fprintf(fh," %E\t%E\t%E\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
d7639a	464
SP	465	fclose(fh);
	466	return TS_SUCCESS;
	467	}
	468
	469
	470	ts_bool fprint_vertex_list(FILE fh,ts_vertex_list vlist){
	471	ts_uint i,j;
	472	for(i=0;i<vlist->n;i++){
8f6a69	473	fprintf(fh," %.17E\t%.17E\t%.17E\t%u\n",vlist->vtx[i]->x,
SP	474	vlist->vtx[i]->y, vlist->vtx[i]->z,
	475	vlist->vtx[i]->neigh_no);
	476	for(j=0;j<vlist->vtx[i]->neigh_no;j++){
	477	fprintf(fh,"\t%u",(ts_uint)(vlist->vtx[i]->neigh[j]->idx));
a6b1b5	478	//-vlist->vtx+1));
d7639a	479	}
SP	480	fprintf(fh,"\n");
	481	}
	482	return TS_SUCCESS;
	483	}
	484
	485	ts_bool fprint_tristar(FILE fh, ts_vesicle vesicle){
	486	ts_uint i,j;
a6b1b5	487	for(i=0;i<vesicle->vlist->n;i++){
8f6a69	488	fprintf(fh,"\t%u",vesicle->vlist->vtx[i]->tristar_no);
SP	489	for(j=0;j<vesicle->vlist->vtx[i]->tristar_no;j++){
	490	fprintf(fh,"\t%u",(ts_uint)(vesicle->vlist->vtx[i]->tristar[j]->idx));//-vesicle->tlist->tria+1));
d7639a	491	}
SP	492	fprintf(fh,"\n");
	493	}
	494	return TS_SUCCESS;
	495	}
	496
	497	ts_bool fprint_triangle_list(FILE fh, ts_vesicle vesicle){
a6b1b5	498	ts_triangle_list *tlist=vesicle->tlist;
d7639a	499	ts_uint i,j;
SP	500	for(i=0;i<tlist->n;i++){
41a035	501	fprintf(fh,"\t%u",tlist->tria[i]->neigh_no);
SP	502	for(j=0;j<tlist->tria[i]->neigh_no;j++){
	503	fprintf(fh,"\t%u",(ts_uint)(tlist->tria[i]->neigh[j]->idx));//-tlist->tria+1));
d7639a	504	}
SP	505	fprintf(fh,"\n");
	506	for(j=0;j<3;j++){
41a035	507	fprintf(fh,"\t%u",(ts_uint)(tlist->tria[i]->vertex[j]->idx));//-vesicle->vlist->vtx+1));
d7639a	508	}
SP	509	fprintf(fh,"\n");
41a035	510	fprintf(fh,"%.17E\t%.17E\t%.17E\n",tlist->tria[i]->xnorm,
SP	511	tlist->tria[i]->ynorm,tlist->tria[i]->znorm);
d7639a	512	fprintf(fh,"0.00000000000000000\n0.00000000000000000\n");
SP	513	}
	514	return TS_SUCCESS;
	515	}
	516
	517	ts_bool fprint_vertex_data(FILE fh,ts_vertex_list vlist){
	518	ts_uint i,j;
	519	for(i=0;i<vlist->n;i++){
	520	fprintf(fh," %.17E\t%.17E\t%.17E\t%.17E\t%.17E\t%u\n",
8f6a69	521	vlist->vtx[i]->xk,vlist->vtx[i]->c,vlist->vtx[i]->energy,
SP	522	vlist->vtx[i]->energy_h, vlist->vtx[i]->curvature, 0);
	523	for(j=0;j<vlist->vtx[i]->neigh_no;j++){
	524	fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length_dual);
d7639a	525	}
SP	526	fprintf(fh,"\n");
8f6a69	527	for(j=0;j<vlist->vtx[i]->neigh_no;j++){
SP	528	fprintf(fh," %.17E", vlist->vtx[i]->bond[j]->bond_length);
d7639a	529	}
SP	530	fprintf(fh,"\n");
	531	}
	532	return TS_SUCCESS;
	533	}
	534
	535	ts_bool fprint_bonds(FILE fh,ts_vesicle vesicle){
	536	ts_uint i;
a6b1b5	537	for(i=0;i<vesicle->blist->n;i++){
e016c4	538	fprintf(fh,"\t%u\t%u\n",(ts_uint)(vesicle->blist->bond[i]->vtx1->idx),
a6b1b5	539	//-vesicle->vlist->vtx+1),
e016c4	540	(ts_uint)(vesicle->blist->bond[i]->vtx2->idx));
a6b1b5	541	//-vesicle->vlist.vtx+1));
d7639a	542	}
SP	543	return TS_SUCCESS;
	544	}
	545
	546
	547	ts_bool write_dout_fcompat_file(ts_vesicle vesicle, ts_char filename){
	548	FILE *fh;
	549	fh=fopen(filename, "w");
	550	if(fh==NULL){
	551	err("Cannot open file %s for writing");
	552	return TS_FAIL;
	553	}
	554	fprintf(fh,"%.17E\n%.17E\n",vesicle->stepsize,vesicle->dmax);
a6b1b5	555	fprint_vertex_list(fh,vesicle->vlist);
d7639a	556	fprint_tristar(fh,vesicle);
SP	557	fprint_triangle_list(fh,vesicle);
a6b1b5	558	fprint_vertex_data(fh,vesicle->vlist);
d7639a	559	fprint_bonds(fh,vesicle);
SP	560	fclose(fh);
	561	return TS_SUCCESS;
	562	}
	563
	564	ts_bool read_tape_fcompat_file(ts_vesicle vesicle, ts_char filename){
	565	FILE *fh;
	566	char line[255];
	567	fh=fopen(filename, "r");
	568	if(fh==NULL){
	569	err("Cannot open file for reading... Nonexistant file?");
	570	return TS_FAIL;
	571	}
a6b1b5	572	ts_uint retval=1;
d7639a	573	while(retval!=EOF){
SP	574	retval=fscanf(fh,"%s",line);
	575
	576	fprintf(stderr,"%s",line);
	577	}
	578	fclose(fh);
	579	return TS_SUCCESS;
	580	}
	581
	582	ts_bool write_master_xml_file(ts_char *filename){
	583	FILE *fh;
	584	ts_char i,j;
	585	ts_uint tstep;
	586	ts_char *number;
	587	fh=fopen(filename, "w");
	588	if(fh==NULL){
	589	err("Cannot open file %s for writing");
	590	return TS_FAIL;
	591	}
	592
	593	fprintf(fh,"<?xml version=\"1.0\"?>\n<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n<Collection>");
	594	DIR *dir = opendir(".");
	595	if(dir){
	596	struct dirent *ent;
	597	tstep=0;
	598	while((ent = readdir(dir)) != NULL)
	599	{
	600	i=rindex(ent->d_name,'.');
	601	if(i==NULL) continue;
	602	if(strcmp(i+1,"vtu")==0){
	603	j=rindex(ent->d_name,'_');
	604	if(j==NULL) continue;
	605	number=strndup(j+1,j-i);
a6b1b5	606	fprintf(fh,"<DataSet timestep=\"%u\" group=\"\" part=\"0\" file=\"%s\"/>\n",atoi(number),ent->d_name);
d7639a	607	tstep++;
SP	608	free(number);
	609	}
	610	}
	611	}
f74313	612	free(dir);
d7639a	613	fprintf(fh,"</Collection>\n</VTKFile>\n");
SP	614	fclose(fh);
	615	return TS_SUCCESS;
	616	}
	617
	618	ts_bool write_vertex_xml_file(ts_vesicle *vesicle, ts_uint timestepno){
a6b1b5	619	ts_vertex_list *vlist=vesicle->vlist;
SP	620	ts_bond_list *blist=vesicle->blist;
	621	ts_vertex **vtx=vlist->vtx;
40aa5b	622	ts_uint i,j;
d7639a	623	char filename[255];
SP	624	FILE *fh;
	625
	626	sprintf(filename,"timestep_%.6u.vtu",timestepno);
	627	fh=fopen(filename, "w");
	628	if(fh==NULL){
	629	err("Cannot open file %s for writing");
	630	return TS_FAIL;
	631	}
	632	/* Here comes header of the file */
40aa5b	633
SP	634	//find number of extra vtxs and bonds of polymeres
3c1ac1	635	ts_uint monono=0, polyno=0, poly_idx=0;
40aa5b	636	ts_bool poly=0;
SP	637	if(vesicle->poly_list!=NULL){
	638	if(vesicle->poly_list->poly[0]!=NULL){
	639	polyno=vesicle->poly_list->n;
	640	monono=vesicle->poly_list->poly[0]->vlist->n;
	641	poly=1;
	642	}
	643	}
d7639a	644	fprintf(fh, "<?xml version=\"1.0\"?>\n<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n <UnstructuredGrid>\n");
40aa5b	645	fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+mononopolyno, blist->n+mononopolyno);
d7639a	646	fprintf(fh,"<PointData Scalars=\"scalars\">\n<DataArray type=\"Int64\" Name=\"scalars\" format=\"ascii\">");
SP	647	for(i=0;i<vlist->n;i++){
a6b1b5	648	fprintf(fh,"%u ",vtx[i]->idx);
d7639a	649	}
40aa5b	650	//polymeres
SP	651	if(poly){
3c1ac1	652	poly_idx=vlist->n;
40aa5b	653	for(i=0;i<vesicle->poly_list->n;i++){
3c1ac1	654	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){
SP	655	fprintf(fh,"%u ", poly_idx);
40aa5b	656	}
SP	657	}
	658	}
d7639a	659
SP	660	fprintf(fh,"</DataArray>\n</PointData>\n<CellData>\n</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n");
	661	for(i=0;i<vlist->n;i++){
8f6a69	662	fprintf(fh,"%e %e %e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
40aa5b	663	}
SP	664	//polymeres
	665	if(poly){
	666	for(i=0;i<vesicle->poly_list->n;i++){
	667	for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
	668	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 );
	669	}
	670	}
d7639a	671	}
SP	672
	673	fprintf(fh,"</DataArray>\n</Points>\n<Cells>\n<DataArray type=\"Int64\" Name=\"connectivity\" format=\"ascii\">");
	674	for(i=0;i<blist->n;i++){
e016c4	675	fprintf(fh,"%u %u\n",blist->bond[i]->vtx1->idx,blist->bond[i]->vtx2->idx);
d7639a	676	}
40aa5b	677	//polymeres
SP	678	if(poly){
3c1ac1	679	poly_idx=vlist->n;
40aa5b	680	for(i=0;i<vesicle->poly_list->n;i++){
SP	681	for(j=0;j<vesicle->poly_list->poly[i]->blist->n;j++){
3c1ac1	682	// 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	683	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	684	}
SP	685	//grafted bonds
3c1ac1	686	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	687	}
SP	688
	689	}
	690
	691
d7639a	692	fprintf(fh,"</DataArray>\n<DataArray type=\"Int64\" Name=\"offsets\" format=\"ascii\">");
40aa5b	693	for (i=2;i<(blist->n+mononopolyno)2+1;i+=2){
d7639a	694	fprintf(fh,"%u ",i);
SP	695	}
	696	fprintf(fh,"\n");
	697	fprintf(fh,"</DataArray>\n<DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">\n");
40aa5b	698	for (i=0;i<blist->n+monono*polyno;i++){
d7639a	699	fprintf(fh,"3 ");
SP	700	}
	701
	702	fprintf(fh,"</DataArray>\n</Cells>\n</Piece>\n</UnstructuredGrid>\n</VTKFile>\n");
	703	fclose(fh);
	704	return TS_SUCCESS;
	705
	706	}
	707
	708	ts_bool write_vertex_vtk_file(ts_vesicle vesicle,ts_char filename, ts_char *text){
a6b1b5	709	ts_vertex_list *vlist=vesicle->vlist;
SP	710	ts_bond_list *blist=vesicle->blist;
	711	ts_vertex **vtx=vlist->vtx;
	712	ts_uint i;
d7639a	713	FILE *fh;
SP	714
	715	fh=fopen(filename, "w");
	716	if(fh==NULL){
	717	err("Cannot open file %s for writing");
	718	return TS_FAIL;
	719	}
	720	/* Here comes header of the file */
	721	// fprintf(stderr,"NSHELL=%u\n",nshell);
	722	fprintf(fh, "# vtk DataFile Version 2.0\n");
	723	/* TODO: Do a sanity check on text. Max 255 char, must not me \n terminated */
	724	fprintf(fh, "%s\n", text);
	725	fprintf(fh,"ASCII\n");
	726	fprintf(fh,"DATASET UNSTRUCTURED_GRID\n");
	727	fprintf(fh,"POINTS %u double\n", vlist->n);
	728	for(i=0;i<vlist->n;i++){
8f6a69	729	fprintf(fh,"%e %e %e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z);
d7639a	730	}
SP	731
	732	fprintf(fh,"CELLS %u %u\n",blist->n,3*blist->n);
	733	for(i=0;i<blist->n;i++){
e016c4	734	fprintf(fh,"2 %u %u\n",blist->bond[i]->vtx1->idx,blist->bond[i]->vtx2->idx);
d7639a	735	}
SP	736	fprintf(fh,"CELL_TYPES %u\n",blist->n);
	737	for(i=0;i<blist->n;i++)
	738	fprintf(fh,"3\n");
	739
	740	fprintf(fh,"POINT_DATA %u\n", vlist->n);
	741	fprintf(fh,"SCALARS scalars long 1\n");
	742	fprintf(fh,"LOOKUP_TABLE default\n");
	743
	744	for(i=0;i<vlist->n;i++)
8f6a69	745	fprintf(fh,"%u\n",vtx[i]->idx);
d7639a	746
SP	747	fclose(fh);
	748	return TS_SUCCESS;
	749	}
	750
	751
	752
d7a113	753	ts_vesicle parsetape(ts_uint mcsweeps, ts_uint inititer, ts_uint iterations){
a2db52	754	long int nshell=17,ncxmax=60, ncymax=60, nczmax=60, npoly=10, nmono=20; // THIS IS DUE TO CONFUSE BUG!
b14a8d	755	char buf=malloc(255sizeof(char));
SP	756	long int brezveze0=1;
	757	long int brezveze1=1;
	758	long int brezveze2=1;
48bb92	759	ts_double xk0=25.0, dmax=1.67,stepsize=0.15,kspring=800.0;
d7a113	760	long int iter=1000, init=1000, mcsw=1000;
40aa5b	761
SP	762
d7639a	763	cfg_opt_t opts[] = {
SP	764	CFG_SIMPLE_INT("nshell", &nshell),
a2db52	765	CFG_SIMPLE_INT("npoly", &npoly),
SP	766	CFG_SIMPLE_INT("nmono", &nmono),
d7639a	767	CFG_SIMPLE_FLOAT("dmax", &dmax),
SP	768	CFG_SIMPLE_FLOAT("xk0",&xk0),
48bb92	769	CFG_SIMPLE_FLOAT("k_spring",&kspring),
d7639a	770	CFG_SIMPLE_FLOAT("stepsize",&stepsize),
SP	771	CFG_SIMPLE_INT("nxmax", &ncxmax),
	772	CFG_SIMPLE_INT("nymax", &ncymax),
	773	CFG_SIMPLE_INT("nzmax", &nczmax),
d7a113	774	CFG_SIMPLE_INT("iterations",&iter),
SP	775	CFG_SIMPLE_INT("mcsweeps",&mcsw),
	776	CFG_SIMPLE_INT("inititer", &init),
d7639a	777	CFG_SIMPLE_BOOL("quiet",&quiet),
314f2d	778	CFG_SIMPLE_STR("multiprocessing",buf),
b14a8d	779	CFG_SIMPLE_INT("smp_cores",&brezveze0),
SP	780	CFG_SIMPLE_INT("cluster_nodes",&brezveze1),
	781	CFG_SIMPLE_INT("distributed_processes",&brezveze2),
d7639a	782	CFG_END()
SP	783	};
	784	cfg_t *cfg;
	785	ts_uint retval;
	786	cfg = cfg_init(opts, 0);
a6b1b5	787	retval=cfg_parse(cfg, "tape");
d7639a	788	if(retval==CFG_FILE_ERROR){
a6b1b5	789	fatal("No tape file.",100);
d7639a	790	}
SP	791	else if(retval==CFG_PARSE_ERROR){
	792	fatal("Invalid tape!",100);
	793	}
b14a8d	794	ts_vesicle *vesicle;
d7a113	795	*iterations=iter;
SP	796	*inititer=init;
	797	*mcsweeps=mcsw;
b14a8d	798	vesicle=initial_distribution_dipyramid(nshell,ncxmax,ncymax,nczmax,stepsize);
a2db52	799	vesicle->poly_list=init_poly_list(npoly,nmono, vesicle->vlist);
48bb92	800	vesicle->spring_constant=kspring;
M	801	poly_assign_spring_const(vesicle);
	802
a2db52	803
d7639a	804	vesicle->nshell=nshell;
SP	805	vesicle->dmax=dmax*dmax;
	806	vesicle->bending_rigidity=xk0;
	807	vesicle->stepsize=stepsize;
a6b1b5	808	vesicle->clist->ncmax[0]=ncxmax;
SP	809	vesicle->clist->ncmax[1]=ncymax;
	810	vesicle->clist->ncmax[2]=nczmax;
	811	vesicle->clist->max_occupancy=8;
b14a8d	812
d7639a	813	cfg_free(cfg);
b14a8d	814	free(buf);
SP	815	// fprintf(stderr,"NSHELL=%u\n",vesicle->nshell);
40aa5b	816
SP	817
	818
b14a8d	819	return vesicle;
d7639a	820
SP	821	}
f74313	822
SP	823
	824	ts_bool read_geometry_file(char fname, ts_vesicle vesicle){
	825	FILE *fh;
	826	ts_uint i, nvtx,nedges,ntria;
	827	ts_uint vtxi1,vtxi2;
	828	float x,y,z;
	829	ts_vertex_list *vlist;
	830	fh=fopen(fname, "r");
	831	if(fh==NULL){
	832	err("Cannot open file for reading... Nonexistant file?");
	833	return TS_FAIL;
	834	}
	835	ts_uint retval;
	836	retval=fscanf(fh,"%u %u %u",&nvtx, &nedges, &ntria);
	837	vesicle->vlist=init_vertex_list(nvtx);
	838	vlist=vesicle->vlist;
	839	for(i=0;i<nvtx;i++){
8f6a69	840	// fscanf(fh,"%F %F %F",&vlist->vtx[i]->x,&vlist->vtx[i]->y,&vlist->vtx[i]->z);
f74313	841	retval=fscanf(fh,"%F %F %F",&x,&y,&z);
8f6a69	842	vlist->vtx[i]->x=x;
SP	843	vlist->vtx[i]->y=y;
	844	vlist->vtx[i]->z=z;
f74313	845	}
SP	846	for(i=0;i<nedges;i++){
	847	retval=fscanf(fh,"%u %u",&vtxi1,&vtxi2);
	848	bond_add(vesicle->blist,vesicle->vlist->vtx[vtxi1-1],vesicle->vlist->vtx[vtxi2-1]);
	849	}
	850	//TODO: neighbours from bonds,
	851	//TODO: triangles from neigbours
	852
	853	// Don't need to read triangles. Already have enough data
	854	/*
	855	for(i=0;i<ntria;i++){
	856	retval=fscanf(fh,"%u %u %u", &bi1, &bi2, &bi3);
8f6a69	857	vtxi1=vesicle->blist->vertex1->idx;
SP	858	vtxi2=vesicle->blist->vertex1->idx;
f74313	859
SP	860	}
	861	*/
41a035	862	if(retval);
f74313	863	fclose(fh);
SP	864
	865
	866
	867	return TS_SUCCESS;
	868	}