From f977e8fe32b6a1a2ae85662ab35fd7c9f647cc3e Mon Sep 17 00:00:00 2001 From: Samo Penic <samo.penic@gmail.com> Date: Mon, 02 Jul 2018 10:50:58 +0000 Subject: [PATCH] Added a funciton that calculates everything needed for spherical harmonics --- src/io.c | 116 ++++++++++++++++++++++++++++++++++++++++++++++++++++----- 1 files changed, 105 insertions(+), 11 deletions(-) diff --git a/src/io.c b/src/io.c index 39dc1e0..3ea3662 100644 --- a/src/io.c +++ b/src/io.c @@ -1,3 +1,4 @@ + /* vim: set ts=4 sts=4 sw=4 noet : */ #include "general.h" #include<stdio.h> @@ -423,7 +424,7 @@ vesicle->tape=parsetape(command_line_args.tape_fullfilename); // recreating space for cells // vesicle->clist=init_cell_list(vesicle->tape->ncxmax, vesicle->tape->ncymax, vesicle->tape->nczmax, vesicle->tape->stepsize); - vesicle->clist->max_occupancy=8; + vesicle->clist->max_occupancy=16; // vesicle->tape=(ts_tape *)malloc(sizeof(ts_tape)); // retval=fread(vesicle->tape, sizeof(ts_tape),1,fh); retval=fread(iteration,sizeof(ts_uint),1,fh); @@ -820,11 +821,12 @@ return TS_SUCCESS; } -ts_bool write_vertex_xml_file(ts_vesicle *vesicle, ts_uint timestepno){ +ts_bool write_vertex_xml_file(ts_vesicle *vesicle, ts_uint timestepno, ts_cluster_list *cstlist){ ts_vertex_list *vlist=vesicle->vlist; ts_bond_list *blist=vesicle->blist; ts_vertex **vtx=vlist->vtx; ts_uint i,j; + //ts_double senergy=0.0; char filename[10000]; char just_name[255]; FILE *fh; @@ -861,8 +863,8 @@ fprintf(fh, "<?xml version=\"1.0\"?>\n<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\" compressor=\"vtkZLibDataCompressor\">\n"); xml_trisurf_data(fh,vesicle); fprintf(fh, " <UnstructuredGrid>\n"); - fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+monono*polyno+fonono*filno, blist->n+monono*polyno+filno*(fonono-1)); - fprintf(fh,"<PointData Scalars=\"scalars\">\n<DataArray type=\"Int64\" Name=\"scalars\" format=\"ascii\">"); + fprintf(fh, "<Piece NumberOfPoints=\"%u\" NumberOfCells=\"%u\">\n",vlist->n+monono*polyno+fonono*filno, blist->n+monono*polyno+filno*(fonono-1)+vesicle->tlist->n); + fprintf(fh,"<PointData Scalars=\"vertices_idx\">\n<DataArray type=\"Int64\" Name=\"vertices_idx\" format=\"ascii\">"); for(i=0;i<vlist->n;i++){ fprintf(fh,"%u ",vtx[i]->idx); } @@ -887,7 +889,40 @@ } fprintf(fh,"</DataArray>\n"); - + if(cstlist!=NULL){ + fprintf(fh,"<DataArray type=\"Int64\" Name=\"vertices_in_cluster\" format=\"ascii\">"); + for(i=0;i<vlist->n;i++){ + if(vtx[i]->cluster!=NULL){ + fprintf(fh,"%u ",vtx[i]->cluster->nvtx); + } else { + fprintf(fh,"-1 "); + } + } + //polymeres + if(poly){ + poly_idx=vlist->n; + for(i=0;i<vesicle->poly_list->n;i++){ + for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){ + fprintf(fh,"-1 "); + } + } + } + //filaments + if(fil){ + poly_idx=vlist->n+monono*polyno; + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){ + // fprintf(stderr,"was here\n"); + fprintf(fh,"-1 "); + } + } + } + + fprintf(fh,"</DataArray>\n"); + + + } + //here comes additional data as needed. Currently only spontaneous curvature fprintf(fh,"<DataArray type=\"Float64\" Name=\"spontaneous_curvature\" format=\"ascii\">"); for(i=0;i<vlist->n;i++){ @@ -914,9 +949,52 @@ } fprintf(fh,"</DataArray>\n"); + //here comes additional data. Energy! + fprintf(fh,"<DataArray type=\"Float64\" Name=\"bending_energy\" format=\"ascii\">"); + for(i=0;i<vlist->n;i++){ + fprintf(fh,"%.17e ",vtx[i]->energy*vtx[i]->xk); + } + //polymeres + if(poly){ + poly_idx=vlist->n; + for(i=0;i<vesicle->poly_list->n;i++){ + for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++,poly_idx++){ + fprintf(fh,"%.17e ", vesicle->poly_list->poly[i]->vlist->vtx[j]->energy* vesicle->poly_list->poly[i]->k); + } + } + } + //filaments + if(fil){ + poly_idx=vlist->n+monono*polyno; + for(i=0;i<vesicle->filament_list->n;i++){ + for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++,poly_idx++){ + // fprintf(stderr,"was here\n"); + fprintf(fh,"%.17e ", vesicle->filament_list->poly[i]->vlist->vtx[j]->energy* vesicle->filament_list->poly[i]->k); + } + } + } + fprintf(fh,"</DataArray>\n"); - fprintf(fh,"</PointData>\n<CellData>\n</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n"); + fprintf(fh,"</PointData>\n<CellData>\n"); + + if(vesicle->tape->stretchswitch==1){ + fprintf(fh,"<DataArray type=\"Float64\" Name=\"stretching_energy\" format=\"ascii\">"); + for(i=0;i<blist->n;i++){ + fprintf(fh, "0.0 "); + } + for(i=0;i<monono*polyno+filno*(fonono-1);i++){ + fprintf(fh,"0.0 "); + } + for(i=0;i<vesicle->tlist->n;i++){ + fprintf(fh,"%.17e ",vesicle->tlist->tria[i]->energy); + } + fprintf(fh,"</DataArray>\n"); + } + + + + fprintf(fh,"</CellData>\n<Points>\n<DataArray type=\"Float64\" Name=\"Koordinate tock\" NumberOfComponents=\"3\" format=\"ascii\">\n"); for(i=0;i<vlist->n;i++){ fprintf(fh,"%.17e %.17e %.17e\n",vtx[i]->x,vtx[i]->y, vtx[i]->z); } @@ -967,17 +1045,24 @@ } } - + for(i=0;i<vesicle->tlist->n;i++){ + fprintf(fh,"%u %u %u\n", vesicle->tlist->tria[i]->vertex[0]->idx, vesicle->tlist->tria[i]->vertex[1]->idx, vesicle->tlist->tria[i]->vertex[2]->idx); + } fprintf(fh,"</DataArray>\n<DataArray type=\"Int64\" Name=\"offsets\" format=\"ascii\">"); for (i=2;i<(blist->n+monono*polyno+(fonono-1)*filno)*2+1;i+=2){ fprintf(fh,"%u ",i); } + for(j=i+1;j<i+3*(vesicle->tlist->n);j+=3){ //let's continue counting from where we left of + fprintf(fh,"%u ", j); + } fprintf(fh,"\n"); fprintf(fh,"</DataArray>\n<DataArray type=\"UInt8\" Name=\"types\" format=\"ascii\">\n"); - for (i=0;i<blist->n+monono*polyno+fonono*filno;i++){ + for (i=0;i<blist->n+monono*polyno+(fonono-1)*filno;i++){ fprintf(fh,"3 "); } - + for(i=0;i<vesicle->tlist->n;i++){ + fprintf(fh,"5 "); + } fprintf(fh,"</DataArray>\n</Cells>\n</Piece>\n</UnstructuredGrid>\n</VTKFile>\n"); fclose(fh); return TS_SUCCESS; @@ -1086,6 +1171,7 @@ CFG_SIMPLE_INT("nmono", &tape->nmono), CFG_SIMPLE_INT("nfil",&tape->nfil), CFG_SIMPLE_INT("nfono",&tape->nfono), + CFG_SIMPLE_INT("internal_poly",&tape->internal_poly), CFG_SIMPLE_INT("R_nucleus",&tape->R_nucleus), CFG_SIMPLE_FLOAT("R_nucleusX",&tape->R_nucleusX), CFG_SIMPLE_FLOAT("R_nucleusY",&tape->R_nucleusY), @@ -1097,6 +1183,8 @@ CFG_SIMPLE_INT("constvolswitch",&tape->constvolswitch), CFG_SIMPLE_INT("constareaswitch",&tape->constareaswitch), CFG_SIMPLE_FLOAT("constvolprecision",&tape->constvolprecision), + CFG_SIMPLE_INT("stretchswitch",&tape->stretchswitch), + CFG_SIMPLE_FLOAT("xkA0",&tape->xkA0), CFG_SIMPLE_FLOAT("pressure",&tape->pressure), CFG_SIMPLE_FLOAT("k_spring",&tape->kspring), CFG_SIMPLE_FLOAT("xi",&tape->xi), @@ -1107,12 +1195,16 @@ CFG_SIMPLE_INT("iterations",&tape->iterations), CFG_SIMPLE_INT("mcsweeps",&tape->mcsweeps), CFG_SIMPLE_INT("inititer", &tape->inititer), - CFG_SIMPLE_BOOL("quiet",&tape->quiet), - CFG_SIMPLE_STR("multiprocessing",tape->multiprocessing), + CFG_SIMPLE_BOOL("quiet",(cfg_bool_t *)&tape->quiet), + CFG_SIMPLE_STR("multiprocessing",&tape->multiprocessing), CFG_SIMPLE_INT("smp_cores",&tape->brezveze0), CFG_SIMPLE_INT("cluster_nodes",&tape->brezveze1), CFG_SIMPLE_INT("distributed_processes",&tape->brezveze2), CFG_SIMPLE_INT("spherical_harmonics_coefficients",&tape->shc), + CFG_SIMPLE_INT("number_of_vertices_with_c0", &tape->number_of_vertices_with_c0), + CFG_SIMPLE_FLOAT("c0",&tape->c0), + CFG_SIMPLE_FLOAT("w",&tape->w), + CFG_SIMPLE_FLOAT("F",&tape->F), CFG_END() }; cfg_t *cfg; @@ -1147,6 +1239,8 @@ ts_bool getcmdline_tape(cfg_t *cfg, char *opts){ char *commands, *backup, *saveptr, *saveopptr, *command, *operator[2]; + operator[0]=0; + operator[1]=0; ts_uint i,j; commands=(char *)malloc(10000*sizeof(char)); backup=commands; //since the pointer to commands will be lost, we acquire a pointer that will serve as backup. -- Gitblit v1.9.3