| | |
| | | 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 |
| | |
| | | #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; |
| | |
| | | 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; |
| | |
| | | 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); |
| | |
| | | #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 |
| | |
| | | #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; |
| | |
| | | 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; |
| | | } |
| | |
| | | #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; |
| | |
| | | #ifndef _H_FRAME |
| | | #define _H_FRAME |
| | | ts_bool centermass(ts_vesicle *vesicle); |
| | | ts_bool cell_occupation(ts_vesicle *vesicle); |
| | | |
| | | #endif |
| | |
| | | #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); |
| | |
| | | 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; |
| | | } |
| | |
| | | #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> |
| | |
| | | } |
| | | } |
| | | } |
| | | free(dir); |
| | | fprintf(fh,"</Collection>\n</VTKFile>\n"); |
| | | fclose(fh); |
| | | return TS_SUCCESS; |
| | |
| | | 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; |
| | | } |
| | |
| | | |
| | | #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??? |
| | |
| | | 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; |
| | | } |
| | |
| | | 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 |