Added stretching variables

Samo Penic

2018-12-11 6c274defd904dedd6cb89c002f5eeced6c079d47

 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;
}