From e984829db39b2778e4f66c34524329ad09749c45 Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Mon, 11 Jul 2016 19:29:21 +0000 Subject: [PATCH] Added possibility of internal pegs. It can break the system however --- src/dumpstate.c | 112 ++++++++++++++++++++++++++++++++++++++++++++++++++++--- 1 files changed, 105 insertions(+), 7 deletions(-) diff --git a/src/dumpstate.c b/src/dumpstate.c index 4d7dde2..4493f71 100644 --- a/src/dumpstate.c +++ b/src/dumpstate.c @@ -1,3 +1,4 @@ +/* vim: set ts=4 sts=4 sw=4 noet : */ #include <string.h> #include "general.h" #include <stdio.h> @@ -6,6 +7,8 @@ #include "dumpstate.h" #include <libxml/parser.h> #include <libxml/tree.h> +#include "vertex.h" +#include "energy.h" ts_vesicle *vtk2vesicle(char *filename, ts_tape *tape){ @@ -46,16 +49,23 @@ cur_node=cur_node->children; while(cur_node!=NULL){ - fprintf(stderr,"Node name is: %s\n",cur_node->name); + //fprintf(stderr,"Node name is: %s\n",cur_node->name); cur_node=cur_node->next; if(strcmp((char *)cur_node->name,"PointData")==0) vtk_index2vesicle(cur_node->children->next->children, vesicle); - if(strcmp((char *)cur_node->name,"Points")==0) break; - if(strcmp((char *)cur_node->name,"Cells")==0) break; - + if(strcmp((char *)cur_node->name,"Points")==0) vtk_coordinates(cur_node->children->next->children,vesicle); + if(strcmp((char *)cur_node->name,"Cells")==0) { + vtk_neighbours(cur_node->children->next->children,vesicle); + break; //segfaults, because it finds another cells. Why? + } } - - - + //we have vertices and neighbour relations all set, but unordered. Let's do bonds, ordering, triangles... + vesicle->vlist = vtk_sort_neighbours(vesicle->blist,vesicle->vlist); + init_triangles(vesicle); + init_triangle_neighbours(vesicle); + init_common_vertex_triangle_neighbours(vesicle); + init_normal_vectors(vesicle->tlist); + mean_curvature_and_energy(vesicle); + ts_fprintf(stdout,"restoring from vtk dump finished!\n"); return TS_SUCCESS; } @@ -72,3 +82,91 @@ //fprintf(stderr,"idx[11]=%d\n",vesicle->vlist->vtx[11]->idx); return TS_SUCCESS; } + +ts_bool vtk_coordinates(xmlNode *node, ts_vesicle *vesicle){ + //fprintf(stderr, "vsebina: %s\n",node->content); + ts_uint i; + char *token; + token = strtok((char *)node->content, " "); + for(i=0;i<vesicle->vlist->n;i++){ + vesicle->vlist->vtx[i]->x=atof(token); +// fprintf(stderr,"%e, ",atof(token)); + token=strtok(NULL," "); + vesicle->vlist->vtx[i]->y=atof(token); + token=strtok(NULL,"\n"); + vesicle->vlist->vtx[i]->z=atof(token); + if(i>0) token=strtok(NULL," "); + } + return TS_SUCCESS; + +} + +ts_bool vtk_neighbours(xmlNode *node, ts_vesicle *vesicle){ + // fprintf(stderr, "vsebina: %s\n",node->content); + ts_uint i; + ts_uint idx1, idx2; + char *token; + token = strtok((char *)node->content, " "); + for(i=0;i<3*(vesicle->vlist->n-2);i++){ + idx1=atoi(token); + token=strtok(NULL,"\n"); + idx2=atoi(token); + vtx_add_neighbour(vesicle->vlist->vtx[idx1], vesicle->vlist->vtx[idx2]); + vtx_add_neighbour(vesicle->vlist->vtx[idx2], vesicle->vlist->vtx[idx1]); + fprintf(stderr, "%d. povezujem %d in %d\n",i,idx1,idx2); + if(i>0) token=strtok(NULL," "); + } + return TS_SUCCESS; +} + +/* this is almost exact copy of init_sort_neighbours */ +ts_vertex_list *vtk_sort_neighbours(ts_bond_list *blist,ts_vertex_list *vlist){ + ts_vertex **vtx=vlist->vtx -1; // take a look at dipyramid function for comment. + ts_uint i,l,j,jj,jjj,k=0; + ts_double eps=0.001; // Take a look if EPS from math.h can be used + +/*lets initialize memory for temporary vertex_list. Should we write a function instead */ + ts_vertex_list *tvlist=vertex_list_copy(vlist); + ts_vertex **tvtx=tvlist->vtx -1; /* again to compensate for 0-indexing */ + + ts_double dist2; // Square of distance of neighbours + ts_double direct; // Something, dont know what, but could be normal of some kind + for(i=1;i<=vlist->n;i++){ + k++; // WHY i IS NOT GOOD?? + vtx_add_cneighbour(blist,tvtx[k], tvtx[vtx[i]->neigh[0]->idx+1]); //always add 1st + jjj=1; + jj=1; + for(l=2;l<=vtx[i]->neigh_no;l++){ + for(j=2;j<=vtx[i]->neigh_no;j++){ + dist2=vtx_distance_sq(vtx[i]->neigh[j-1],vtx[i]->neigh[jj-1]); + direct=vtx_direct(vtx[i],vtx[i]->neigh[j-1],vtx[i]->neigh[jj-1]); +// TODO: check if fabs can be used with all floating point types!! + if( (direct>0.0) && (j!=jjj) ){ + vtx_add_cneighbour(blist,tvtx[k],tvtx[vtx[i]->neigh[j-1]->idx+1]); + jjj=jj; + jj=j; + // break; + } + } + } + if(vtx[i]->neigh_no!=tvtx[i]->neigh_no){ + fprintf(stderr,"Hej! pri vtx=%d se razikuje stevilo sosedov (prej=%d, potem=%d)\n",i, vtx[i]->neigh_no, tvtx[i]->neigh_no); + } + } + if(eps>dist2); +/* We use the temporary vertex for our main vertices and we abandon main + * vertices, because their neighbours are not correctly ordered */ + // tvtx=vlist->vtx; + // vlist->vtx=tvtx; + // tvlist->vtx=vtx; + vtx_list_free(vlist); +/* Let's make a check if the number of bonds is correct */ + if((blist->n)!=3*(tvlist->n-2)){ + ts_fprintf(stderr,"Number of bonds is %u should be %u!\n", blist->n, 3*(tvlist->n-2)); + fatal("Number of bonds is not 3*(no_vertex-2).",4); + } + + return tvlist; +} + + -- Gitblit v1.9.3