From f74313919463dd08d93faeefd40900e2917c0157 Mon Sep 17 00:00:00 2001
From: Samo Penic <samo@andromeda>
Date: Sun, 05 Dec 2010 09:00:47 +0000
Subject: [PATCH] Some rewritting done.
---
src/Makefile.am | 2
src/frame.h | 6 +
src/vertex.c | 34 ++++++++
src/io.c | 53 ++++++++++++
src/tape | 11 ++
src/energy.c | 81 ++++++++++----------
src/initial_distribution.c | 2
src/frame.c | 36 ++++----
src/vertex.h | 5 +
9 files changed, 169 insertions(+), 61 deletions(-)
diff --git a/src/Makefile.am b/src/Makefile.am
index 024f250..e19ba39 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -1,5 +1,5 @@
trisurfdir=../
trisurf_PROGRAMS=trisurf
-trisurf_SOURCES=general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c main.c
+trisurf_SOURCES=general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c frame.c energy.c main.c
trisurf_LDFLAGS = -lm -lconfuse
trisurf_CFLAGS = -Wall -Werror
diff --git a/src/energy.c b/src/energy.c
index 61720ee..540686b 100644
--- a/src/energy.c
+++ b/src/energy.c
@@ -6,14 +6,13 @@
#include<stdio.h>
ts_bool mean_curvature_and_energy(ts_vesicle *vesicle){
- ts_uint i, jj, jjp;
+ ts_uint i;
- ts_vertex_list *vlist=&vesicle->vlist;
- ts_vertex *vtx=vlist->vertex;
+ ts_vertex_list *vlist=vesicle->vlist;
+ ts_vertex **vtx=vlist->vtx;
for(i=0;i<vlist->n;i++){
- //should call with zero index!!!
- energy_vertex(&vtx[i]);
+ energy_vertex(vtx[i]);
}
return TS_SUCCESS;
@@ -23,6 +22,7 @@
inline ts_bool energy_vertex(ts_vertex *vtx){
// ts_vertex *vtx=&vlist->vertex[n]-1; // Caution! 0 Indexed value!
// ts_triangle *tristar=vtx->tristar-1;
+ ts_vertex_data *data=vtx->data;
ts_uint jj;
ts_uint jjp,jjm;
ts_vertex *j,*jp, *jm;
@@ -30,20 +30,20 @@
ts_double s=0,xh=0,yh=0,zh=0,txn=0,tyn=0,tzn=0;
ts_double x1,x2,x3,ctp,ctm,tot,xlen;
ts_double h,ht;
- for(jj=1; jj<=vtx->neigh_no;jj++){
+ for(jj=1; jj<=data->neigh_no;jj++){
jjp=jj+1;
- if(jjp>vtx->neigh_no) jjp=1;
+ if(jjp>data->neigh_no) jjp=1;
jjm=jj-1;
- if(jjm<1) jjm=vtx->neigh_no;
- j=vtx->neigh[jj-1];
- jp=vtx->neigh[jjp-1];
- jm=vtx->neigh[jjm-1];
- jt=vtx->tristar[jj-1];
- x1=vertex_distance_sq(vtx,jp); //shouldn't be zero!
- x2=vertex_distance_sq(j,jp); // shouldn't be zero!
- x3=(j->x-jp->x)*(vtx->x-jp->x)+
- (j->y-jp->y)*(vtx->y-jp->y)+
- (j->z-jp->z)*(vtx->z-jp->z);
+ if(jjm<1) jjm=data->neigh_no;
+ j=data->neigh[jj-1];
+ jp=data->neigh[jjp-1];
+ jm=data->neigh[jjm-1];
+ jt=data->tristar[jj-1];
+ x1=vtx_distance_sq(vtx,jp); //shouldn't be zero!
+ x2=vtx_distance_sq(j,jp); // shouldn't be zero!
+ x3=(j->data->x-jp->data->x)*(data->x-jp->data->x)+
+ (j->data->y-jp->data->y)*(data->y-jp->data->y)+
+ (j->data->z-jp->data->z)*(data->z-jp->data->z);
#ifdef TS_DOUBLE_DOUBLE
ctp=x3/sqrt(x1*x2-x3*x3);
@@ -54,11 +54,11 @@
#ifdef TS_DOUBLE_LONGDOUBLE
ctp=x3/sqrtl(x1*x2-x3*x3);
#endif
- x1=vertex_distance_sq(vtx,jm);
- x2=vertex_distance_sq(j,jm);
- x3=(j->x-jm->x)*(vtx->x-jm->x)+
- (j->y-jm->y)*(vtx->y-jm->y)+
- (j->z-jm->z)*(vtx->z-jm->z);
+ x1=vtx_distance_sq(vtx,jm);
+ x2=vtx_distance_sq(j,jm);
+ x3=(j->data->x-jm->data->x)*(data->x-jm->data->x)+
+ (j->data->y-jm->data->y)*(data->y-jm->data->y)+
+ (j->data->z-jm->data->z)*(data->z-jm->data->z);
#ifdef TS_DOUBLE_DOUBLE
ctm=x3/sqrt(x1*x2-x3*x3);
#endif
@@ -70,26 +70,26 @@
#endif
tot=ctp+ctm;
tot=0.5*tot;
- xlen=vertex_distance_sq(j,vtx);
+ xlen=vtx_distance_sq(j,vtx);
#ifdef TS_DOUBLE_DOUBLE
- vtx->bond_length[jj-1]=sqrt(xlen);
+ data->bond[jj-1]->data->bond_length=sqrt(xlen);
#endif
#ifdef TS_DOUBLE_FLOAT
- vtx->bond_length[jj-1]=sqrtf(xlen);
+ data->bond[jj-1]->data->bond_length=sqrtf(xlen);
#endif
#ifdef TS_DOUBLE_LONGDOUBLE
- vtx->bond_length[jj-1]=sqrtl(xlen);
+ data->bond[jj-1]->data->bond_length=sqrtl(xlen);
#endif
- vtx->bond_length_dual[jj-1]=tot*vtx->bond_length[jj-1];
+ data->bond[jj-1]->data->bond_length_dual=tot*data->bond[jj-1]->data->bond_length;
s+=tot*xlen;
- xh+=tot*(j->x - vtx->x);
- yh+=tot*(j->y - vtx->y);
- zh+=tot*(j->z - vtx->z);
- txn+=jt->xnorm;
- tyn+=jt->ynorm;
- tzn+=jt->znorm;
+ xh+=tot*(j->data->x - data->x);
+ yh+=tot*(j->data->y - data->y);
+ zh+=tot*(j->data->z - data->z);
+ txn+=jt->data->xnorm;
+ tyn+=jt->data->ynorm;
+ tzn+=jt->data->znorm;
}
h=xh*xh+yh*yh+zh*zh;
@@ -97,26 +97,27 @@
s=s/4.0;
#ifdef TS_DOUBLE_DOUBLE
if(ht>=0.0) {
- vtx->curvature=sqrt(h);
+ data->curvature=sqrt(h);
} else {
- vtx->curvature=-sqrt(h);
+ data->curvature=-sqrt(h);
}
#endif
#ifdef TS_DOUBLE_FLOAT
if(ht>=0.0) {
- vtx->curvature=sqrtf(h);
+ data->curvature=sqrtf(h);
} else {
- vtx->curvature=-sqrtf(h);
+ data->curvature=-sqrtf(h);
}
#endif
#ifdef TS_DOUBLE_LONGDOUBLE
if(ht>=0.0) {
- vtx->curvature=sqrtl(h);
+ data->curvature=sqrtl(h);
} else {
- vtx->curvature=-sqrtl(h);
+ data->curvature=-sqrtl(h);
}
#endif
- vtx->energy=0.5*s*(vtx->curvature/s-vtx->c)*(vtx->curvature/s-vtx->c);
+//TODO: MAJOR!!!! What is vtx->data->c?????????????? Here it is 0!
+ data->energy=0.5*s*(data->curvature/s-data->c)*(data->curvature/s-data->c);
return TS_SUCCESS;
}
diff --git a/src/frame.c b/src/frame.c
index f89f6bc..a578188 100644
--- a/src/frame.c
+++ b/src/frame.c
@@ -1,38 +1,38 @@
#include<stdlib.h>
#include "general.h"
#include "cell.h"
-
+#include "frame.h"
ts_bool centermass(ts_vesicle *vesicle){
- ts_uint i;
+ ts_uint i, n=vesicle->vlist->n;
+ ts_vertex **vtx=vesicle->vlist->vtx;
vesicle->cm[0]=0;
vesicle->cm[1]=0;
vesicle->cm[2]=0;
- for(i=0;i<vesicle->vlist.n;i++){
- vesicle->cm[0]+=vesicle->vlist.vertex[i].x;
- vesicle->cm[1]+=vesicle->vlist.vertex[i].y;
- vesicle->cm[2]+=vesicle->vlist.vertex[i].z;
+ for(i=0;i<n;i++){
+ vesicle->cm[0]+=vtx[i]->data->x;
+ vesicle->cm[1]+=vtx[i]->data->y;
+ vesicle->cm[2]+=vtx[i]->data->z;
}
- vesicle->cm[0]/=(float)vesicle->vlist.n;
- vesicle->cm[1]/=(float)vesicle->vlist.n;
- vesicle->cm[2]/=(float)vesicle->vlist.n;
+ vesicle->cm[0]/=(float)n;
+ vesicle->cm[1]/=(float)n;
+ vesicle->cm[2]/=(float)n;
return TS_SUCCESS;
}
-ts_bool cell_ocupation(ts_vesicle *vesicle){
- ts_uint i,j,cellidx;
+ts_bool cell_occupation(ts_vesicle *vesicle){
+ ts_uint i,cellidx, n=vesicle->vlist->n;
ts_double shift;
ts_double dcell;
- shift=(ts_double) vesicle->clist.ncmax[0]/2;
+ shift=(ts_double) vesicle->clist->ncmax[0]/2;
dcell=1.0/(1.0 + vesicle->stepsize);
- ts_uint ncx, ncy,ncz;
- cell_list_cell_ocupation_clear(&vesicle->clist);
- for(i=0;i<vesicle->vlist.n;i++){
+ cell_list_cell_occupation_clear(vesicle->clist);
+ for(i=0;i<n;i++){
- cellidx=vertex_self_avoidance(vesicle, &vesicle->vlist.vertex[i]);
- vesicle->vlist.vertex[i].cell=&(vesicle->clist.cell[cellidx]);
+ cellidx=vertex_self_avoidance(vesicle, vesicle->vlist->vtx[i]);
+ vesicle->vlist->vtx[i]->data->cell=vesicle->clist->cell[cellidx];
- cell_add_vertex(&vesicle->clist.cell[cellidx],&vesicle->vlist.vertex[i]);
+ cell_add_vertex(vesicle->clist->cell[cellidx],vesicle->vlist->vtx[i]);
// if(ncx > vesicle->clist.ncmax[0]) vesicle->clist.ncmax[0]=ncx;
// if(ncy > vesicle->clist.ncmax[1]) vesicle->clist.ncmax[1]=ncy;
diff --git a/src/frame.h b/src/frame.h
index e69de29..b0e451c 100644
--- a/src/frame.h
+++ b/src/frame.h
@@ -0,0 +1,6 @@
+#ifndef _H_FRAME
+#define _H_FRAME
+ts_bool centermass(ts_vesicle *vesicle);
+ts_bool cell_occupation(ts_vesicle *vesicle);
+
+#endif
diff --git a/src/initial_distribution.c b/src/initial_distribution.c
index d5f1e94..4881e11 100644
--- a/src/initial_distribution.c
+++ b/src/initial_distribution.c
@@ -8,6 +8,7 @@
#include "vertex.h"
#include "triangle.h"
#include "initial_distribution.h"
+#include "energy.h"
ts_vesicle *initial_distribution_dipyramid(ts_uint nshell, ts_uint ncmax1, ts_uint ncmax2, ts_uint ncmax3, ts_double stepsize){
ts_fprintf(stderr,"Starting initial_distribution on vesicle with %u shells!...\n",nshell);
@@ -24,6 +25,7 @@
retval = init_triangles(vesicle);
retval = init_triangle_neighbours(vesicle);
retval = init_common_vertex_triangle_neighbours(vesicle);
+ retval = mean_curvature_and_energy(vesicle);
ts_fprintf(stderr,"initial_distribution finished!\n");
return vesicle;
}
diff --git a/src/io.c b/src/io.c
index 86d5d02..1a7a527 100644
--- a/src/io.c
+++ b/src/io.c
@@ -1,7 +1,9 @@
-#include<general.h>
+#include "general.h"
#include<stdio.h>
#include "io.h"
-#include "confuse.h"
+#include <confuse.h>
+#include "vertex.h"
+#include "bond.h"
#include<string.h>
#include<stdlib.h>
#include <sys/types.h>
@@ -194,6 +196,7 @@
}
}
}
+ free(dir);
fprintf(fh,"</Collection>\n</VTKFile>\n");
fclose(fh);
return TS_SUCCESS;
@@ -330,3 +333,49 @@
return TS_SUCCESS;
}
+
+
+ts_bool read_geometry_file(char *fname, ts_vesicle *vesicle){
+ FILE *fh;
+ ts_uint i, nvtx,nedges,ntria;
+ ts_uint vtxi1,vtxi2;
+ float x,y,z;
+ ts_vertex_list *vlist;
+ fh=fopen(fname, "r");
+ if(fh==NULL){
+ err("Cannot open file for reading... Nonexistant file?");
+ return TS_FAIL;
+ }
+ ts_uint retval;
+ retval=fscanf(fh,"%u %u %u",&nvtx, &nedges, &ntria);
+ vesicle->vlist=init_vertex_list(nvtx);
+ vlist=vesicle->vlist;
+ for(i=0;i<nvtx;i++){
+ // fscanf(fh,"%F %F %F",&vlist->vtx[i]->data->x,&vlist->vtx[i]->data->y,&vlist->vtx[i]->data->z);
+ retval=fscanf(fh,"%F %F %F",&x,&y,&z);
+ vlist->vtx[i]->data->x=x;
+ vlist->vtx[i]->data->y=y;
+ vlist->vtx[i]->data->z=z;
+ }
+ for(i=0;i<nedges;i++){
+ retval=fscanf(fh,"%u %u",&vtxi1,&vtxi2);
+ bond_add(vesicle->blist,vesicle->vlist->vtx[vtxi1-1],vesicle->vlist->vtx[vtxi2-1]);
+ }
+ //TODO: neighbours from bonds,
+ //TODO: triangles from neigbours
+
+// Don't need to read triangles. Already have enough data
+ /*
+ for(i=0;i<ntria;i++){
+ retval=fscanf(fh,"%u %u %u", &bi1, &bi2, &bi3);
+ vtxi1=vesicle->blist->data->vertex1->idx;
+ vtxi2=vesicle->blist->data->vertex1->idx;
+
+ }
+ */
+ fclose(fh);
+
+
+
+ return TS_SUCCESS;
+}
diff --git a/src/tape b/src/tape
index 26b4d44..64f529a 100644
--- a/src/tape
+++ b/src/tape
@@ -21,3 +21,14 @@
#shut up if we are using cluster!!!
quiet=false
+
+#what type of multiprocessing? (*none, smp, cluster, distributed, cuda, auto)
+#currently only none makes sense.
+multiprocessing=none
+#how many cores are allowed to process in SMP?
+smp_cores=2
+#how many nodes in cluster?
+cluster_nodes=50
+#max number of processes in distributed (voluntary) environment
+distributed_processes=50
+#cuda options???
diff --git a/src/vertex.c b/src/vertex.c
index 379d43e..dfb0155 100644
--- a/src/vertex.c
+++ b/src/vertex.c
@@ -214,3 +214,37 @@
return TS_SUCCESS;
}
+
+/* ****************************************************************** */
+/* ***** New vertex copy operations. Inherently they are slow. ***** */
+/* ****************************************************************** */
+
+ts_vertex *vtx_copy(ts_vertex *ovtx){
+ ts_vertex *cvtx=(ts_vertex *)malloc(sizeof(ts_vertex));
+ memcpy((void *)cvtx,(void *)ovtx,sizeof(ts_vertex));
+ cvtx->data=(ts_vertex_data *)malloc(sizeof(ts_vertex_data));
+ memcpy((void *)cvtx->data,(void *)ovtx->data,sizeof(ts_vertex_data));
+ cvtx->data->neigh=NULL;
+ cvtx->data->neigh_no=0;
+ cvtx->data->tristar_no=0;
+ cvtx->data->bond_no=0;
+ cvtx->data->tristar=NULL;
+ cvtx->data->bond=NULL;
+ cvtx->data->cell=NULL;
+ return cvtx;
+}
+//TODO: needs to be done
+ts_vertex **vtx_neigh_copy(ts_vertex_list *vlist,ts_vertex *ovtx){
+ return NULL;
+}
+
+ts_vertex_list *vertex_list_copy(ts_vertex_list *ovlist){
+ ts_uint i;
+ ts_vertex_list *vlist=(ts_vertex_list *)malloc(sizeof(ts_vertex_list));
+ memcpy((void *)vlist, (void *)ovlist, sizeof(ts_vertex_list));
+ ts_vertex **vtx=(ts_vertex **)malloc(vlist->n*sizeof(ts_vertex *));
+ for(i=0;i<vlist->n;i++){
+ vtx[i]=vtx_copy(ovlist->vtx[i]);
+ }
+ return vlist;
+}
diff --git a/src/vertex.h b/src/vertex.h
index 4e22dd2..0b353aa 100644
--- a/src/vertex.h
+++ b/src/vertex.h
@@ -26,4 +26,9 @@
ts_bool vtx_set_global_values(ts_vesicle *vesicle);
inline ts_double vtx_direct(ts_vertex *vtx1, ts_vertex *vtx2, ts_vertex *vtx3);
inline ts_bool vertex_add_tristar(ts_vertex *vtx, ts_triangle *tristarmem);
+
+ts_vertex *vtx_copy(ts_vertex *ovtx);
+ts_vertex **vtx_neigh_copy(ts_vertex_list *vlist,ts_vertex *ovtx);
+ts_vertex_list *vertex_list_copy(ts_vertex_list *ovlist);
+
#endif
--
Gitblit v1.9.3