From ecdd71929b38a504e2b0ae1a9c105a5e1ac28e4a Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@fe.uni-lj.si> Date: Mon, 29 Feb 2016 08:56:06 +0000 Subject: [PATCH] Dump and restore seems equal to original in paraview. --- src/sh.c | 52 +++++++++++++++++++++++++++++++++++++++++++++++++--- 1 files changed, 49 insertions(+), 3 deletions(-) diff --git a/src/sh.c b/src/sh.c index 0c4b06d..9cb623e 100644 --- a/src/sh.c +++ b/src/sh.c @@ -2,7 +2,8 @@ #include<stdlib.h> #include "general.h" #include "sh.h" - +#include "io.h" +#include <string.h> ts_spharm *sph_init(ts_vertex_list *vlist, ts_uint l){ @@ -49,6 +50,7 @@ ts_bool 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->sumUlm2[i]!=NULL) free(sph->sumUlm2[i]); @@ -189,7 +191,7 @@ K=-sqrt(1.0/(M_PI))*cos(m*fi); } - return K*sqrt((2.0*l+1.0)/2.0*fac1/fac2)*plgndr(l,abs(m),cos(theta)); + return K*sqrt((2.0*l+1.0)/2.0*(ts_double)(fac1/fac2))*plgndr(l,abs(m),cos(theta)); } @@ -200,7 +202,7 @@ #ifdef TS_DOUBLE_DOUBLE coord->e1=sqrt(x*x+y*y+z*z); if(z==0) coord->e3=M_PI/2.0; - else coord->e3=atan(sqrt(x*x+y*y)/z); + else coord->e3=atan2(sqrt(x*x+y*y),z); coord->e2=atan2(y,x); #endif #ifdef TS_DOUBLE_FLOAT @@ -218,6 +220,23 @@ return TS_SUCCESS; } + + +ts_bool sph2cart(ts_coord *coord){ + coord->coord_type=TS_COORD_CARTESIAN; + ts_double x,y,z; + + x=coord->e1*cos(coord->e2)*sin(coord->e3); + y=coord->e1*sin(coord->e2)*sin(coord->e3); + z=coord->e1*cos(coord->e3); + + coord->e1=x; + coord->e2=y; + coord->e3=z; + + return TS_SUCCESS; +} + /* Function returns radius of the sphere with the same volume as vesicle (r0) */ ts_double getR0(ts_vesicle *vesicle){ @@ -379,3 +398,30 @@ sph->N++; return TS_SUCCESS; } + + +ts_bool saveAvgUlm2(ts_vesicle *vesicle){ + + FILE *fh; + char filename[10000]; + strcpy(filename, command_line_args.path); + strcat(filename, "sph2out.dat"); + fh=fopen(filename, "w"); + if(fh==NULL){ + err("Cannot open file %s for writing"); + return TS_FAIL; + } + + ts_spharm *sph=vesicle->sphHarmonics; + ts_int i,j; + fprintf(fh,"l,\tm,\tulm^2avg\n"); + for(i=0;i<sph->l;i++){ + for(j=0;j<2*i+1;j++){ + fprintf(fh,"%d,\t%d,\t%e\n", i, j-i, sph->sumUlm2[i][j]/(ts_double)sph->N); + + } + fprintf(fh,"\n"); + } + fclose(fh); + return TS_SUCCESS; +} -- Gitblit v1.9.3