From 62b11daa07f943c5ab047049cf2e2803188d467c Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@fe.uni-lj.si> Date: Tue, 05 Jul 2016 10:11:56 +0000 Subject: [PATCH] Added stub function for clustering --- src/io.c | 581 ++++++++++++++++++++++++++++++++++++++++++++++++++------- 1 files changed, 503 insertions(+), 78 deletions(-) diff --git a/src/io.c b/src/io.c index d4a9cc4..7066657 100644 --- a/src/io.c +++ b/src/io.c @@ -1,7 +1,7 @@ +/* vim: set ts=4 sts=4 sw=4 noet : */ #include "general.h" #include<stdio.h> #include "io.h" -#include <confuse.h> #include "vertex.h" #include "bond.h" #include<string.h> @@ -10,22 +10,23 @@ #include <dirent.h> #include "initial_distribution.h" #include "poly.h" - +#include "cell.h" #include <getopt.h> #include <sys/stat.h> #include <sys/types.h> #include <dirent.h> #include <errno.h> +#include <snapshot.h> /** DUMP STATE TO DISK DRIVE **/ ts_bool dump_state(ts_vesicle *vesicle, ts_uint iteration){ /* save current state with wrong pointers. Will fix that later */ ts_uint i,j,k; - FILE *fh=fopen("dump.bin","wb"); + FILE *fh=fopen(command_line_args.dump_fullfilename,"wb"); /* dump vesicle */ - fwrite(vesicle, sizeof(ts_vesicle),1,fh); + fwrite(vesicle, sizeof(ts_vesicle)-sizeof(ts_double),1,fh); /* dump vertex list */ fwrite(vesicle->vlist, sizeof(ts_vertex_list),1,fh); /* dump bond list */ @@ -36,6 +37,8 @@ fwrite(vesicle->clist, sizeof(ts_cell_list),1,fh); /* dump poly list */ fwrite(vesicle->poly_list, sizeof(ts_poly_list),1,fh); + /* dump filament list */ + fwrite(vesicle->filament_list, sizeof(ts_poly_list),1,fh); /* level 1 complete */ /*dump vertices*/ @@ -118,6 +121,43 @@ } } + + /*dump filamentes grandes svinjas */ + for(i=0;i<vesicle->filament_list->n;i++){ + fwrite(vesicle->filament_list->poly[i],sizeof(ts_poly),1,fh); + fwrite(vesicle->filament_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh); + fwrite(vesicle->filament_list->poly[i]->blist,sizeof(ts_bond_list),1,fh); + } + + /* dump filamentes vertex(monomer) list*/ + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){ + fwrite(vesicle->filament_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh); + /* dump offset for neigh and bond */ + for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh_no;k++){ + // off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->neigh[k]-vesicle->poly_list->poly[i]->vlist->vtx[0]); + fwrite(&vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh[k]->idx,sizeof(ts_uint),1,fh); + } + for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->bond_no;k++){ + //off=(ts_ulong)(vesicle->poly_list->poly[i]->vlist->vtx[j]->bond[k]-vesicle->poly_list->poly[i]->blist->bond[0]); + fwrite(&vesicle->filament_list->poly[i]->vlist->vtx[j]->bond[k]->idx,sizeof(ts_uint),1,fh); + } + } + } + /* dump poly bonds between monomers list*/ + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){ + fwrite(vesicle->filament_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh); + /* dump vtx1 and vtx2 offsets */ + //off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx1-vesicle->poly_list->poly[i]->vlist->vtx[0]); + fwrite(&vesicle->filament_list->poly[i]->blist->bond[j]->vtx1->idx,sizeof(ts_uint),1,fh); +// off=(ts_ulong)(vesicle->poly_list->poly[i]->blist->bond[j]->vtx2-vesicle->poly_list->poly[i]->vlist->vtx[0]); + fwrite(&vesicle->filament_list->poly[i]->blist->bond[j]->vtx2->idx,sizeof(ts_uint),1,fh); + } + } + + + /* pointer offsets for fixing the restored pointers */ /* need pointers for vlist->vtx @@ -130,8 +170,9 @@ poly_list->poly->bond */ - fwrite(vesicle->clist, sizeof(ts_cell_list),1, fh); - +// fwrite(vesicle->clist, sizeof(ts_cell_list),1, fh); +/* write tape information on vesicle */ +// fwrite(vesicle->tape,sizeof(ts_tape),1,fh); fwrite(&iteration, sizeof(ts_uint),1,fh); fclose(fh); return TS_SUCCESS; @@ -141,14 +182,28 @@ /** RESTORE DUMP FROM DISK **/ ts_vesicle *restore_state(ts_uint *iteration){ ts_uint i,j,k; - FILE *fh=fopen("dump.bin","rb"); + FILE *fh=fopen(command_line_args.dump_fullfilename,"rb"); + + struct stat sb; + if (stat(command_line_args.dump_fullfilename, &sb) == -1) { + //dump file does not exist. + return NULL; + } + + //check if it is regular file + if((sb.st_mode & S_IFMT) != S_IFREG) { + //dump file is not a regular file. + ts_fprintf(stderr,"Dump file is not a regular file!\n"); + return NULL; + } + ts_uint retval; ts_uint idx; /* we restore all the data from the dump */ /* restore vesicle */ ts_vesicle *vesicle=(ts_vesicle *)calloc(1,sizeof(ts_vesicle)); - retval=fread(vesicle, sizeof(ts_vesicle),1,fh); + retval=fread(vesicle, sizeof(ts_vesicle)-sizeof(ts_double),1,fh); // fprintf(stderr,"was here! %e\n",vesicle->dmax); /* restore vertex list */ @@ -166,6 +221,9 @@ /* restore poly list */ vesicle->poly_list=(ts_poly_list *)calloc(1,sizeof(ts_poly_list)); retval=fread(vesicle->poly_list, sizeof(ts_poly_list),1,fh); + /* restore filament list */ + vesicle->filament_list=(ts_poly_list *)calloc(1,sizeof(ts_poly_list)); + retval=fread(vesicle->filament_list, sizeof(ts_poly_list),1,fh); /* level 1 complete */ /* prerequisity. Bonds must be malloced before vertexes are recreated */ @@ -308,15 +366,66 @@ } } -// recreating space for cells // - vesicle->clist=(ts_cell_list *)malloc(sizeof(ts_cell_list)); - retval=fread(vesicle->clist, sizeof(ts_cell_list), 1,fh); - vesicle->clist->cell=(ts_cell **)malloc(sizeof(ts_cell *)*vesicle->clist->ncmax[0]*vesicle->clist->ncmax[1]*vesicle->clist->ncmax[2]); - for(i=0;i<vesicle->clist->ncmax[0]*vesicle->clist->ncmax[1]*vesicle->clist->ncmax[2];i++){ - vesicle->clist->cell[i]=(ts_cell *)calloc(1,sizeof(ts_cell)); - vesicle->clist->cell[i]->idx=i+1; // We enumerate cells! Probably never required! - } + /*restore filaments */ + vesicle->filament_list->poly = (ts_poly **)calloc(vesicle->filament_list->n,sizeof(ts_poly *)); + for(i=0;i<vesicle->filament_list->n;i++){ + vesicle->filament_list->poly[i]=(ts_poly *)calloc(1,sizeof(ts_poly)); + retval=fread(vesicle->filament_list->poly[i],sizeof(ts_poly),1,fh); + vesicle->filament_list->poly[i]->vlist=(ts_vertex_list *)calloc(1,sizeof(ts_vertex_list)); + retval=fread(vesicle->filament_list->poly[i]->vlist,sizeof(ts_vertex_list),1,fh); + vesicle->filament_list->poly[i]->blist=(ts_bond_list *)calloc(1,sizeof(ts_bond_list)); + retval=fread(vesicle->filament_list->poly[i]->blist,sizeof(ts_bond_list),1,fh); + /* initialize adress space for pointers that will hold specific vertices (monomers) and bonds */ + vesicle->filament_list->poly[i]->vlist->vtx=(ts_vertex **)calloc(vesicle->filament_list->poly[i]->vlist->n,sizeof(ts_vertex *)); + vesicle->filament_list->poly[i]->blist->bond=(ts_bond **)calloc(vesicle->filament_list->poly[i]->blist->n,sizeof(ts_bond *)); + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){ + vesicle->filament_list->poly[i]->vlist->vtx[j]=(ts_vertex *)malloc(sizeof(ts_vertex)); + } + for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){ + vesicle->filament_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond)); + } + } + + + /* restore poly vertex(monomer) list*/ + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){ + retval=fread(vesicle->filament_list->poly[i]->vlist->vtx[j],sizeof(ts_vertex),1,fh); + + /* restore neigh and bonds */ + 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 *)); + for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh_no;k++){ + retval=fread(&idx,sizeof(ts_uint),1,fh); + vesicle->filament_list->poly[i]->vlist->vtx[j]->neigh[k]=vesicle->filament_list->poly[i]->vlist->vtx[idx]; + } + 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 *)); + for(k=0;k<vesicle->filament_list->poly[i]->vlist->vtx[j]->bond_no;k++){ + retval=fread(&idx,sizeof(ts_uint),1,fh); + vesicle->filament_list->poly[i]->vlist->vtx[j]->bond[k]=vesicle->filament_list->poly[i]->blist->bond[idx]; + } + + } + } + + /* restore poly bonds between monomers list*/ + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){ + // vesicle->poly_list->poly[i]->blist->bond[j]=(ts_bond *)malloc(sizeof(ts_bond)); + retval=fread(vesicle->filament_list->poly[i]->blist->bond[j],sizeof(ts_bond),1,fh); + /* restore vtx1 and vtx2 */ + retval=fread(&idx,sizeof(ts_uint),1,fh); + vesicle->filament_list->poly[i]->blist->bond[j]->vtx1=vesicle->filament_list->poly[i]->vlist->vtx[idx]; + retval=fread(&idx,sizeof(ts_uint),1,fh); + vesicle->filament_list->poly[i]->blist->bond[j]->vtx2=vesicle->filament_list->poly[i]->vlist->vtx[idx]; + } + } + vesicle->tape=parsetape(command_line_args.tape_fullfilename); +// recreating space for cells // + vesicle->clist=init_cell_list(vesicle->tape->ncxmax, vesicle->tape->ncymax, vesicle->tape->nczmax, vesicle->tape->stepsize); + vesicle->clist->max_occupancy=8; +// vesicle->tape=(ts_tape *)malloc(sizeof(ts_tape)); +// retval=fread(vesicle->tape, sizeof(ts_tape),1,fh); retval=fread(iteration,sizeof(ts_uint),1,fh); if(retval); fclose(fh); @@ -326,9 +435,15 @@ ts_bool parse_args(int argc, char **argv){ - int c, retval; - DIR *dir; - path[0]=0; + int c, retval; + struct stat sb; + sprintf(command_line_args.path, "./"); //clear string; + sprintf(command_line_args.output_fullfilename,"output.pvd"); + sprintf(command_line_args.dump_fullfilename,"dump.bin"); + sprintf(command_line_args.tape_fullfilename,"tape"); + sprintf(command_line_args.tape_templatefull,"./tape"); + FILE *file; + while (1) { static struct option long_options[] = @@ -336,15 +451,19 @@ {"force-from-tape", no_argument, &(command_line_args.force_from_tape), 1}, {"reset-iteration-count", no_argument, &(command_line_args.reset_iteration_count), 1}, {"tape", no_argument, 0, 't'}, + {"version", no_argument, 0, 'v'}, {"output-file", required_argument, 0, 'o'}, {"directory", required_argument, 0, 'd'}, - {"dump-file", required_argument,0, 'f'}, + {"dump-filename", required_argument,0, 'f'}, + {"tape-options",required_argument,0,'c'}, + {"tape-template", required_argument,0,0}, + {"restore-from-vtk",required_argument,0,0}, {0, 0, 0, 0} }; /* getopt_long stores the option index here. */ int option_index = 0; - c = getopt_long (argc, argv, "d:fot", + c = getopt_long (argc, argv, "d:f:o:t:c:v", long_options, &option_index); /* Detect the end of the options. */ @@ -357,60 +476,63 @@ /* If this option set a flag, do nothing else now. */ if (long_options[option_index].flag != 0) break; - printf ("option %s", long_options[option_index].name); +/* printf ("option %s", long_options[option_index].name); if (optarg) - printf (" with arg %s", optarg); - printf ("\n"); + printf (" with arg %s", optarg); + printf ("\n"); */ + //TODO: find a better way. + if(strcmp(long_options[option_index].name,"tape-template")==0){ + strcpy(command_line_args.tape_templatefull,optarg); + } + if(strcmp(long_options[option_index].name,"restore-from-vtk")==0){ + strcpy(command_line_args.dump_from_vtk,optarg); + } + break; + case 'v': + fprintf(stdout,"TRISURF-NG v. %s, compiled on: %s %s.\n", TS_VERSION, __DATE__, __TIME__); + fprintf(stdout,"Programming done by: Samo Penic and Miha Fosnaric\n"); + fprintf(stdout,"Released under terms of GPLv3\n"); + exit(0); + + case 'c': + strcpy(command_line_args.tape_opts,optarg); + break; + case 't': //tape + strcpy(command_line_args.tape_fullfilename,optarg); break; - case 't': - //check if tape exists. If not, fail immediately. - puts ("option -t\n"); - break; - - case 'o': - //set filename of master output file - printf ("option -o with value `%s'\n", optarg); - break; + case 'o': //set filename of master pvd output file + strcpy(command_line_args.output_fullfilename, optarg); + break; case 'd': //check if directory exists. If not create one. If creation is //successful, set directory for output files. //printf ("option -d with value `%s'\n", optarg); - dir = opendir(optarg); - if (dir) - { - /* Directory exists. */ - closedir(dir); - } - else if (ENOENT == errno) - { - /* Directory does not exist. */ + if (stat(optarg, &sb) == -1) { + //directory does not exist retval=mkdir(optarg, 0700); if(retval){ - fatal("Could not create requested directory. Check if you have permissions",1); + fatal("Could not create requested directory. Check if you have permissions",1); } } - else - { - /* opendir() failed for some other reason. */ - fatal("Could not check if directory exists. Reason unknown",1); + //check if is a proper directory + else if((sb.st_mode & S_IFMT) != S_IFDIR) { + //it is not a directory. fatal error. + ts_fprintf(stderr,"%s is not a directory!\n",optarg); + fatal("Cannot continue",1); } - ts_fprintf(stdout,"\n*** Using output directory: %s\n\n", optarg); -// sprintf(path,"%s", optarg); - strcpy(path, optarg); - // ts_fprintf(stdout,"ok!\n"); - + strcpy(command_line_args.path, optarg); break; case 'f': - //check if dump file specified exists. Defaults to dump.bin + strcpy(command_line_args.dump_fullfilename, optarg); break; case '?': /* getopt_long already printed an error message. */ - - ts_fprintf(stderr,"\n\nhere comes the help.\n\n"); + print_help(stdout); +//ts_fprintf(stderr,"\n\nhere comes the help.\n\n"); fatal("Ooops, read help first",1); break; @@ -419,10 +541,90 @@ } } +//Here we set correct values for full filenames! + char *buffer=(char *)malloc(10000*sizeof(char)); + //correct the path and add trailing / + if(command_line_args.path[strlen(command_line_args.path)-1]!='/') strcat(command_line_args.path,"/"); + +/* master pvd output file */ + strcpy(buffer,command_line_args.path); + strcat(buffer,command_line_args.output_fullfilename); + if ((file = fopen(buffer, "w")) == NULL) { + fprintf(stderr,"Could not create output file %s!\n", buffer); + fatal("Please specify correct output file or check permissions of the file",1); + //there is a tape template. make a copy into desired directory + + } else { + fclose(file); + strcpy(command_line_args.output_fullfilename,buffer); + } + +/* tape file */ + strcpy(buffer,command_line_args.path); + strcat(buffer,command_line_args.tape_fullfilename); + if (stat(buffer, &sb) == -1) { + + //tape does not exist. does tape template exist? + if(stat(command_line_args.tape_templatefull, &sb)==-1){ + ts_fprintf(stderr,"Tape '%s' does not exist and no tape template was specified (or does not exist)!\n",buffer); + fatal("Please select correct tape or check permissions of the file",1); + } else { + //tape template found + fatal("Samo did not program template copy yet",1); + } + } else { + strcpy(command_line_args.tape_fullfilename,buffer); + } + + +/* dump file */ + strcpy(buffer,command_line_args.path); + strcat(buffer,command_line_args.dump_fullfilename); + //check if dump file exist first. + if (stat(buffer, &sb) == -1) { + //no dump file. check if we can create one. + if ((file = fopen(buffer, "w")) == NULL) { + fprintf(stderr,"Could not create dump file '%s'!\n",buffer); + fatal("Please specify correct dump file or check permissions of the file",1); + } else { + fclose(file); + //good, file is writeable. delete it for now. + remove(buffer); + } + } + strcpy(command_line_args.dump_fullfilename, buffer); + + + free(buffer); return TS_SUCCESS; } + +ts_bool print_help(FILE *fd){ + fprintf(fd,"TRISURF-NG v. %s, compiled on: %s %s.\n", TS_VERSION, __DATE__, __TIME__); + fprintf(fd,"Programming done by: Samo Penic and Miha Fosnaric\n"); + fprintf(fd,"Released under terms of GPLv3\n\n"); + + 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"); + 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"); + 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"); + fprintf(fd,"Flags:\n\n"); + 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"); + 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"); + fprintf(fd,"--reset-iteration-count\t\t starts simulation from the beginning (using binary dump).\n"); + fprintf(fd,"--tape (or -t)\t\t specifies tape filename. For --force-from-tape and restoring from binary dump. Defaults to 'tape'.\n"); + fprintf(fd,"--version (or -v)\t\t Prints version information.\n"); + fprintf(fd,"--output-file (or -o)\t\t Specifies filename of .PVD file. Defaults to 'output.pvd'\n"); + fprintf(fd,"--dump-filename (or -f)\t\t specifies filename for binary dump&restore. Defaults to 'dump.bin'\n\n\n"); + fprintf(fd,"Examples:\n\n"); + fprintf(fd,"trisurf --force-from-tape\n"); + fprintf(fd,"trisurf --reset-iteration-count\n"); + fprintf(fd,"trisurf --restore-from-vtk filename.vtu\n"); + fprintf(fd,"\n\n"); + + return TS_SUCCESS; +} @@ -594,7 +796,7 @@ } fprintf(fh,"<?xml version=\"1.0\"?>\n<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n<Collection>"); - DIR *dir = opendir("."); + DIR *dir = opendir(command_line_args.path); if(dir){ struct dirent *ent; tstep=0; @@ -623,10 +825,13 @@ ts_bond_list *blist=vesicle->blist; ts_vertex **vtx=vlist->vtx; ts_uint i,j; - char filename[255]; + char filename[10000]; + char just_name[255]; FILE *fh; + strcpy(filename,command_line_args.path); + sprintf(just_name,"timestep_%.6u.vtu",timestepno); + strcat(filename,just_name); - sprintf(filename,"timestep_%.6u.vtu",timestepno); fh=fopen(filename, "w"); if(fh==NULL){ err("Cannot open file %s for writing"); @@ -635,8 +840,8 @@ /* Here comes header of the file */ //find number of extra vtxs and bonds of polymeres - ts_uint monono=0, polyno=0, poly_idx=0; - ts_bool poly=0; + ts_uint monono=0, polyno=0, poly_idx=0, filno=0, fonono=0; + ts_bool poly=0, fil=0; if(vesicle->poly_list!=NULL){ if(vesicle->poly_list->poly[0]!=NULL){ polyno=vesicle->poly_list->n; @@ -644,9 +849,20 @@ poly=1; } } - fprintf(fh, "<?xml version=\"1.0\"?>\n<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n <UnstructuredGrid>\n"); - fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+monono*polyno, blist->n+monono*polyno); - fprintf(fh,"<PointData Scalars=\"scalars\">\n<DataArray type=\"Int64\" Name=\"scalars\" format=\"ascii\">"); + + if(vesicle->filament_list!=NULL){ + if(vesicle->filament_list->poly[0]!=NULL){ + filno=vesicle->filament_list->n; + fonono=vesicle->filament_list->poly[0]->vlist->n; + fil=1; + } + } + + fprintf(fh, "<?xml version=\"1.0\"?>\n<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n"); + xml_trisurf_data(fh,vesicle); + fprintf(fh, " <UnstructuredGrid>\n"); + fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+monono*polyno+fonono*filno, blist->n+monono*polyno+filno*(fonono-1)+vesicle->tlist->n); + fprintf(fh,"<PointData Scalars=\"vertices_idx\">\n<DataArray type=\"Int64\" Name=\"vertices_idx\" format=\"ascii\">"); for(i=0;i<vlist->n;i++){ fprintf(fh,"%u ",vtx[i]->idx); } @@ -659,16 +875,90 @@ } } } + //filaments + if(fil){ + poly_idx=vlist->n+monono*polyno; + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){ + // fprintf(stderr,"was here\n"); + fprintf(fh,"%u ", poly_idx); + } + } + } - fprintf(fh,"</DataArray>\n</PointData>\n<CellData>\n</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n"); + fprintf(fh,"</DataArray>\n"); + + //here comes additional data as needed. Currently only spontaneous curvature + fprintf(fh,"<DataArray type=\"Float64\" Name=\"spontaneous_curvature\" format=\"ascii\">"); for(i=0;i<vlist->n;i++){ - fprintf(fh,"%e %e %e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z); + fprintf(fh,"%.17e ",vtx[i]->c); + } + //polymeres + if(poly){ + poly_idx=vlist->n; + for(i=0;i<vesicle->poly_list->n;i++){ + for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){ + fprintf(fh,"%.17e ", vesicle->poly_list->poly[i]->vlist->vtx[j]->c); + } + } + } + //filaments + if(fil){ + poly_idx=vlist->n+monono*polyno; + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){ + // fprintf(stderr,"was here\n"); + fprintf(fh,"%.17e ", vesicle->filament_list->poly[i]->vlist->vtx[j]->c); + } + } + } + fprintf(fh,"</DataArray>\n"); + + //here comes additional data. Energy! + fprintf(fh,"<DataArray type=\"Float64\" Name=\"bending_energy\" format=\"ascii\">"); + for(i=0;i<vlist->n;i++){ + fprintf(fh,"%.17e ",vtx[i]->energy*vtx[i]->xk); + } + //polymeres + if(poly){ + poly_idx=vlist->n; + for(i=0;i<vesicle->poly_list->n;i++){ + for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){ + fprintf(fh,"%.17e ", vesicle->poly_list->poly[i]->vlist->vtx[j]->energy* vesicle->poly_list->poly[i]->k); + } + } + } + //filaments + if(fil){ + poly_idx=vlist->n+monono*polyno; + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){ + // fprintf(stderr,"was here\n"); + fprintf(fh,"%.17e ", vesicle->filament_list->poly[i]->vlist->vtx[j]->energy* vesicle->filament_list->poly[i]->k); + } + } + } + fprintf(fh,"</DataArray>\n"); + + + + fprintf(fh,"</PointData>\n<CellData>\n</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n"); + for(i=0;i<vlist->n;i++){ + fprintf(fh,"%.17e %.17e %.17e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z); } //polymeres if(poly){ for(i=0;i<vesicle->poly_list->n;i++){ for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){ - 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 ); + 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 ); + } + } + } + //filaments + if(fil){ + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){ + 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 ); } } } @@ -691,17 +981,36 @@ } + //filaments + if(fil){ + poly_idx=vlist->n+monono*polyno; + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->blist->n;j++){ + 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); +// fprintf(stderr,"was here\n"); + + } + } + } + for(i=0;i<vesicle->tlist->n;i++){ + 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); + } fprintf(fh,"</DataArray>\n<DataArray type=\"Int64\" Name=\"offsets\" format=\"ascii\">"); - for (i=2;i<(blist->n+monono*polyno)*2+1;i+=2){ + for (i=2;i<(blist->n+monono*polyno+(fonono-1)*filno)*2+1;i+=2){ fprintf(fh,"%u ",i); } + for(j=i+1;j<i+3*(vesicle->tlist->n);j+=3){ //let's continue counting from where we left of + fprintf(fh,"%u ", j); + } fprintf(fh,"\n"); fprintf(fh,"</DataArray>\n<DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">\n"); - for (i=0;i<blist->n+monono*polyno;i++){ + for (i=0;i<blist->n+monono*polyno+(fonono-1)*filno;i++){ fprintf(fh,"3 "); } - + for(i=0;i<vesicle->tlist->n;i++){ + fprintf(fh,"5 "); + } fprintf(fh,"</DataArray>\n</Cells>\n</Piece>\n</UnstructuredGrid>\n</VTKFile>\n"); fclose(fh); return TS_SUCCESS; @@ -753,26 +1062,77 @@ +ts_bool write_pov_file(ts_vesicle *vesicle, char *filename){ + FILE *fh; + ts_uint i; + + fh=fopen(filename, "w"); + if(fh==NULL){ + err("Cannot open file %s for writing"); + return TS_FAIL; + } + + for(i=0;i<vesicle->tlist->n;i++){ + + fprintf(fh,"\ttriangle {"); + fprintf(fh,"\t<%e,%e,%e> <%e,%e,%e> <%e,%e,%e> }\n", + vesicle->tlist->tria[i]->vertex[0]->x, + vesicle->tlist->tria[i]->vertex[0]->y, + vesicle->tlist->tria[i]->vertex[0]->z, + + vesicle->tlist->tria[i]->vertex[1]->x, + vesicle->tlist->tria[i]->vertex[1]->y, + vesicle->tlist->tria[i]->vertex[1]->z, + + vesicle->tlist->tria[i]->vertex[2]->x, + vesicle->tlist->tria[i]->vertex[2]->y, + vesicle->tlist->tria[i]->vertex[2]->z + ); + } + + fclose(fh); + return TS_SUCCESS; +} + + ts_tape *parsetape(char *filename){ - // long int nshell=17,ncxmax=60, ncymax=60, nczmax=60, npoly=10, nmono=20, pswitch=0; // THIS IS DUE TO CONFUSE BUG! + FILE *fd = fopen (filename, "r"); + long length; + size_t size; + fseek (fd, 0, SEEK_END); + length = ftell (fd); + fseek (fd, 0, SEEK_SET); + size=fread (tapetxt, 1, length, fd); + fclose(fd); + if(size); + ts_tape *tape=parsetapebuffer(tapetxt); + return tape; +} + +ts_tape *parsetapebuffer(char *buffer){ ts_tape *tape=(ts_tape *)calloc(1,sizeof(ts_tape)); tape->multiprocessing=calloc(255,sizeof(char)); - /* long int brezveze0=1; - long int brezveze1=1; - long int brezveze2=1; - ts_double xk0=25.0, dmax=1.67,stepsize=0.15,kspring=800.0,pressure=0.0; - long int iter=1000, init=1000, mcsw=1000; -*/ cfg_opt_t opts[] = { CFG_SIMPLE_INT("nshell", &tape->nshell), CFG_SIMPLE_INT("npoly", &tape->npoly), CFG_SIMPLE_INT("nmono", &tape->nmono), - CFG_SIMPLE_FLOAT("dmax", &tape->dmax), + CFG_SIMPLE_INT("nfil",&tape->nfil), + CFG_SIMPLE_INT("nfono",&tape->nfono), + CFG_SIMPLE_INT("R_nucleus",&tape->R_nucleus), + CFG_SIMPLE_FLOAT("R_nucleusX",&tape->R_nucleusX), + CFG_SIMPLE_FLOAT("R_nucleusY",&tape->R_nucleusY), + CFG_SIMPLE_FLOAT("R_nucleusZ",&tape->R_nucleusZ), + CFG_SIMPLE_FLOAT("dmax", &tape->dmax), + CFG_SIMPLE_FLOAT("dmin_interspecies", &tape->dmin_interspecies), CFG_SIMPLE_FLOAT("xk0",&tape->xk0), CFG_SIMPLE_INT("pswitch",&tape->pswitch), + CFG_SIMPLE_INT("constvolswitch",&tape->constvolswitch), + CFG_SIMPLE_INT("constareaswitch",&tape->constareaswitch), + CFG_SIMPLE_FLOAT("constvolprecision",&tape->constvolprecision), CFG_SIMPLE_FLOAT("pressure",&tape->pressure), CFG_SIMPLE_FLOAT("k_spring",&tape->kspring), + CFG_SIMPLE_FLOAT("xi",&tape->xi), CFG_SIMPLE_FLOAT("stepsize",&tape->stepsize), CFG_SIMPLE_INT("nxmax", &tape->ncxmax), CFG_SIMPLE_INT("nymax", &tape->ncymax), @@ -785,12 +1145,16 @@ CFG_SIMPLE_INT("smp_cores",&tape->brezveze0), CFG_SIMPLE_INT("cluster_nodes",&tape->brezveze1), CFG_SIMPLE_INT("distributed_processes",&tape->brezveze2), + CFG_SIMPLE_INT("spherical_harmonics_coefficients",&tape->shc), + CFG_SIMPLE_INT("number_of_vertices_with_c0", &tape->number_of_vertices_with_c0), + CFG_SIMPLE_FLOAT("c0",&tape->c0), + CFG_SIMPLE_FLOAT("w",&tape->w), CFG_END() }; cfg_t *cfg; ts_uint retval; cfg = cfg_init(opts, 0); - retval=cfg_parse(cfg, filename); + retval=cfg_parse_buf(cfg, buffer); if(retval==CFG_FILE_ERROR){ fatal("No tape file.",100); } @@ -798,6 +1162,8 @@ fatal("Invalid tape!",100); } + /* here we override all values read from tape with values from commandline*/ + getcmdline_tape(cfg,command_line_args.tape_opts); cfg_free(cfg); @@ -813,6 +1179,65 @@ } + +ts_bool getcmdline_tape(cfg_t *cfg, char *opts){ + + char *commands, *backup, *saveptr, *saveopptr, *command, *operator[2]; + ts_uint i,j; + commands=(char *)malloc(10000*sizeof(char)); + backup=commands; //since the pointer to commands will be lost, we acquire a pointer that will serve as backup. + strcpy(commands,opts); + for(i=0; ;i++, commands=NULL){ + //breaks comma separated list of commands into specific commands. + command=strtok_r(commands,",",&saveptr); + if(command==NULL) break; +// fprintf(stdout,"Command %d: %s\n",i,command); + //extracts name of command and value of command into operator[2] array. + for(j=0; j<2;j++,command=NULL){ + operator[j]=strtok_r(command,"=",&saveopptr); + if(operator[j]==NULL) break; +// fprintf(stdout," ---> Operator %d: %s\n",j,operator[j]); + } + //1. check: must have 2 operators. + if(j!=2) fatal("Error. Command line tape options are not formatted properly",1); + + // cfg_setstr(cfg,operator[0],operator[1]); + cmdline_to_tape(cfg,operator[0],operator[1]); + //2. check: must be named properly. + //3. check: must be of right format (integer, double, string, ...) + + } + free(backup); + return TS_SUCCESS; +} + + +ts_bool cmdline_to_tape(cfg_t *cfg, char *key, char *val){ + + cfg_opt_t *cfg_opt=cfg_getopt(cfg,key); + if(cfg_opt==NULL) fatal("Commandline tape option not recognised",1); //return TS_FAIL; + switch (cfg_opt->type){ + case CFGT_INT: + cfg_setint(cfg,key,atol(val)); + break; + case CFGT_FLOAT: + cfg_setfloat(cfg,key,atof(val)); + break; +/* case CFGT_BOOL: + cfg_setbool(cfg,operator[0],operator[1]); + break; */ + case CFGT_STR: + cfg_setstr(cfg,key,val); + break; + default: + break; + + } + return TS_SUCCESS; +} + + + ts_bool read_geometry_file(char *fname, ts_vesicle *vesicle){ FILE *fh; ts_uint i, nvtx,nedges,ntria; -- Gitblit v1.9.3