src/Makefile.am | ●●●●● patch | view | raw | blame | history | |
src/general.h | ●●●●● patch | view | raw | blame | history | |
src/initial_distribution.c | ●●●●● patch | view | raw | blame | history | |
src/main.c | ●●●●● patch | view | raw | blame | history | |
src/shcomplex.c | ●●●●● patch | view | raw | blame | history | |
src/shcomplex.h | ●●●●● patch | view | raw | blame | history | |
src/tape | ●●●●● patch | view | raw | blame | history | |
src/timestep.c | ●●●●● patch | view | raw | blame | history | |
src/vesicle.c | ●●●●● patch | view | raw | blame | history |
src/Makefile.am
@@ -1,17 +1,17 @@ trisurfdir=../ trisurf_PROGRAMS = trisurf trisurf_SOURCES = general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c frame.c energy.c timestep.c vertexmove.c bondflip.c main.c poly.c stats.c sh.c trisurf_SOURCES = general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c frame.c energy.c timestep.c vertexmove.c bondflip.c main.c poly.c stats.c sh.c shcomplex.c #trisurf_LDFLAGS = -lm -lconfuse shdiscoverdir=../ shdiscover_PROGRAMS= shdiscover shdiscover_SOURCES= general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c energy.c sh.c shdiscover.c poly.c stats.c shdiscover_SOURCES= general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c energy.c sh.c shdiscover.c poly.c stats.c shcomplex.c AM_CFLAGS = -Wall -Werror co_testdir=../ co_test_PROGRAMS=co_test co_test_SOURCES= general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c energy.c sh.c co_test.c frame.c bondflip.c poly.c stats.c co_test_SOURCES= general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c energy.c sh.c co_test.c frame.c bondflip.c poly.c stats.c shcomplex.c spherical_trisurfdir=../ spherical_trisurf_PROGRAMS = spherical_trisurf spherical_trisurf_SOURCES = general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c frame.c energy.c timestep.c vertexmove.c spherical_trisurf.c sh.c bondflip.c poly.c stats.c spherical_trisurf_SOURCES = general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c frame.c energy.c timestep.c vertexmove.c spherical_trisurf.c sh.c bondflip.c poly.c stats.c shcomplex.c spherical_trisurf_ffdir=../ spherical_trisurf_ff_PROGRAMS = spherical_trisurf_ff spherical_trisurf_ff_SOURCES = general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c frame.c energy.c timestep.c vertexmove.c spherical_trisurf_ff.c sh.c bondflip.c poly.c stats.c spherical_trisurf_ff_SOURCES = general.c vertex.c bond.c triangle.c cell.c vesicle.c initial_distribution.c io.c frame.c energy.c timestep.c vertexmove.c spherical_trisurf_ff.c sh.c bondflip.c poly.c stats.c shcomplex.c src/general.h
@@ -3,6 +3,7 @@ #include<stdarg.h> #include<stdio.h> #include<gsl/gsl_complex.h> /* @brief This is a header file, defining general constants and structures. * @file header.h * @author Samo Penic @@ -222,6 +223,7 @@ typedef struct { ts_uint l; ts_double **ulm; gsl_complex **ulmComplex; ts_double **sumUlm2; ts_uint N; ts_double **co; src/initial_distribution.c
@@ -12,6 +12,7 @@ #include "poly.h" #include "io.h" #include "sh.h" #include "shcomplex.h" ts_vesicle *initial_distribution_dipyramid(ts_uint nshell, ts_uint ncmax1, ts_uint ncmax2, ts_uint ncmax3, ts_double stepsize){ ts_fprintf(stdout,"Starting initial_distribution on vesicle with %u shells!...\n",nshell); @@ -94,7 +95,7 @@ vesicle->pressure= tape->pressure; vesicle->pswitch=tape->pswitch; if(tape->shc>0){ vesicle->sphHarmonics=sph_init(vesicle->vlist,tape->shc); vesicle->sphHarmonics=complex_sph_init(vesicle->vlist,tape->shc); } else { vesicle->sphHarmonics=NULL; src/main.c
@@ -12,6 +12,7 @@ #include "timestep.h" #include "poly.h" #include "sh.h" #include "shcomplex.h" /** Entrance function to the program * @param argv is a number of parameters used in program call (including the program name @@ -53,7 +54,7 @@ vesicle->clist->dmin_interspecies = tape->dmin_interspecies*tape->dmin_interspecies; /* spherical harmonics */ if(tape->shc>0){ vesicle->sphHarmonics=sph_init(vesicle->vlist,tape->shc); vesicle->sphHarmonics=complex_sph_init(vesicle->vlist,tape->shc); } else { vesicle->sphHarmonics=NULL; src/shcomplex.c
New file @@ -0,0 +1,127 @@ #include<math.h> #include<stdlib.h> #include<gsl/gsl_complex.h> #include<gsl/gsl_complex_math.h> #include<gsl/gsl_sf_legendre.h> #include "general.h" #include "sh.h" #include "shcomplex.h" ts_spharm *complex_sph_init(ts_vertex_list *vlist, ts_uint l){ ts_uint j,i; ts_spharm *sph=(ts_spharm *)malloc(sizeof(ts_spharm)); sph->N=0; /* lets initialize Ylm for each vertex. */ sph->Ylmi=(ts_double ***)calloc(l,sizeof(ts_double **)); for(i=0;i<l;i++){ sph->Ylmi[i]=(ts_double **)calloc(2*i+1,sizeof(ts_double *)); for(j=0;j<(2*i+1);j++){ sph->Ylmi[i][j]=(ts_double *)calloc(vlist->n,sizeof(ts_double)); } } /* lets initialize ulm */ sph->ulm=(ts_double **)calloc(l,sizeof(ts_double *)); sph->ulmComplex=(gsl_complex **)calloc(l,sizeof(gsl_complex *)); for(j=0;j<l;j++){ sph->ulm[j]=(ts_double *)calloc(2*j+1,sizeof(ts_double)); sph->ulmComplex[j]=(gsl_complex *)calloc(2*j+1,sizeof(gsl_complex)); } /* lets initialize sum of Ulm2 */ sph->sumUlm2=(ts_double **)calloc(l,sizeof(ts_double *)); for(j=0;j<l;j++){ sph->sumUlm2[j]=(ts_double *)calloc(2*j+1,sizeof(ts_double)); } /* lets initialize co */ //NOTE: C is has zero based indexing. Code is imported from fortran and to comply with original indexes we actually generate one index more. Also second dimension is 2*j+2 instead of 2*j+2. elements starting with 0 are useles and should be ignored! sph->co=(ts_double **)calloc(l+1,sizeof(ts_double *)); for(j=0;j<=l;j++){ sph->co[j]=(ts_double *)calloc(2*j+2,sizeof(ts_double)); } sph->l=l; /* Calculate coefficients that will remain constant during all the simulation */ precomputeShCoeff(sph); return sph; } ts_bool complex_sph_free(ts_spharm *sph){ int i,j; if(sph==NULL) return TS_FAIL; for(i=0;i<sph->l;i++){ if(sph->ulm[i]!=NULL) free(sph->ulm[i]); if(sph->ulmComplex[i]!=NULL) free(sph->ulmComplex[i]); if(sph->sumUlm2[i]!=NULL) free(sph->sumUlm2[i]); if(sph->co[i]!=NULL) free(sph->co[i]); } if(sph->co[sph->l]!=NULL) free(sph->co[sph->l]); if(sph->co != NULL) free(sph->co); if(sph->ulm !=NULL) free(sph->ulm); if(sph->ulmComplex !=NULL) free(sph->ulmComplex); if(sph->Ylmi!=NULL) { for(i=0;i<sph->l;i++){ if(sph->Ylmi[i]!=NULL){ for(j=0;j<i*2+1;j++){ if(sph->Ylmi[i][j]!=NULL) free (sph->Ylmi[i][j]); } free(sph->Ylmi[i]); } } free(sph->Ylmi); } free(sph); return TS_SUCCESS; } ts_bool calculateUlmComplex(ts_vesicle *vesicle){ ts_int i,j,k,m,l; ts_vertex *cvtx; ts_coord coord; /* set all values to zero */ for(i=0;i<vesicle->sphHarmonics->l;i++){ for(j=0;j<2*i+1;j++) GSL_SET_COMPLEX(&(vesicle->sphHarmonics->ulmComplex[i][j]),0.0,0.0); } for(k=0;k<vesicle->vlist->n; k++){ cvtx=vesicle->vlist->vtx[k]; cart2sph(&coord,cvtx->x,cvtx->y,cvtx->z); for(i=0;i<vesicle->sphHarmonics->l;i++){ for(j=0;j<2*i+1;j++){ m=j-i; l=i; if(m>=0){ // fprintf(stderr, "Racunam za l=%d, m=%d\n", l,m); vesicle->sphHarmonics->ulmComplex[i][j]=gsl_complex_add(vesicle->sphHarmonics->ulmComplex[i][j], gsl_complex_conjugate(gsl_complex_mul_real(gsl_complex_polar(1.0,(ts_double)m*coord.e2),cvtx->solAngle*cvtx->relR*gsl_sf_legendre_sphPlm(l,m,cos(coord.e3)))) ); } else { // fprintf(stderr, "Racunam za l=%d, abs(m=%d)\n", l,m); vesicle->sphHarmonics->ulmComplex[i][j]=gsl_complex_add(vesicle->sphHarmonics->ulmComplex[i][j], gsl_complex_conjugate(gsl_complex_mul_real(gsl_complex_polar(1.0,(ts_double)m*coord.e2),cvtx->solAngle*cvtx->relR*pow(-1,m)*gsl_sf_legendre_sphPlm(l,-m,cos(coord.e3)))) ); } } } } return TS_SUCCESS; } ts_bool storeUlmComplex2(ts_vesicle *vesicle){ ts_spharm *sph=vesicle->sphHarmonics; ts_int i,j; for(i=0;i<sph->l;i++){ for(j=0;j<2*i+1;j++){ sph->sumUlm2[i][j]+=gsl_complex_abs2(sph->ulmComplex[i][j]); } } sph->N++; return TS_SUCCESS; } src/shcomplex.h
New file @@ -0,0 +1,8 @@ #ifndef _H_SH_COMPLEX #define _H_SH_COMPLEX #include "general.h" ts_bool storeUlmComplex2(ts_vesicle *vesicle); ts_spharm *complex_sph_init(ts_vertex_list *vlist, ts_uint l); ts_bool complex_sph_free(ts_spharm *sph); ts_bool calculateUlmComplex(ts_vesicle *vesicle); #endif src/tape
@@ -44,15 +44,15 @@ ####### Program Control ############ #how many MC sweeps between subsequent records of states to disk mcsweeps=10 mcsweeps=250 #how many initial mcsweeps*inititer MC sweeps before recording to disk? inititer=0 #how many records do you want on the disk iteration are there in a run? iterations=10 iterations=1000 ###### Spherical harmonics ########### spherical_harmonics_coefficients=21 spherical_harmonics_coefficients=15 #shut up if we are using cluster!!! quiet=false src/timestep.c
@@ -10,6 +10,7 @@ #include "io.h" #include "stats.h" #include "sh.h" #include "shcomplex.h" #include "vesicle.h" ts_bool run_simulation(ts_vesicle *vesicle, ts_uint mcsweeps, ts_uint inititer, ts_uint iterations, ts_uint start_iteration){ @@ -51,9 +52,9 @@ r0=getR0(vesicle); if(vesicle->sphHarmonics!=NULL){ preparationSh(vesicle,r0); calculateYlmi(vesicle); calculateUlm(vesicle); storeUlm2(vesicle); //calculateYlmi(vesicle); calculateUlmComplex(vesicle); storeUlmComplex2(vesicle); saveAvgUlm2(vesicle); fd1=fopen("state.dat","w"); fprintf(fd1,"%e %e\n",vesicle->volume, getR0(vesicle)); src/vesicle.c
@@ -7,6 +7,7 @@ #include "stdlib.h" #include "poly.h" #include "sh.h" #include "shcomplex.h" ts_vesicle *init_vesicle(ts_uint N, ts_uint ncmax1, ts_uint ncmax2, ts_uint ncmax3, ts_double stepsize){ @@ -36,7 +37,7 @@ triangle_list_free(vesicle->tlist); cell_list_free(vesicle->clist); poly_list_free(vesicle->poly_list); sph_free(vesicle->sphHarmonics); complex_sph_free(vesicle->sphHarmonics); free(vesicle); return TS_SUCCESS; }