Trisurf Monte Carlo simulator
Samo Penic
2014-04-06 01d87137bd3f6459388e03a7b076c6f567022030
commit | author | age
d7639a 1 #include<stdlib.h>
SP 2 #include<stdio.h>
aec47d 3 #include<math.h>
SP 4 //#include "io.h"
5 #include "general.h"
6 #include "timestep.h"
7 #include "vertexmove.h"
30ee9c 8 #include "bondflip.h"
d7a113 9 #include "frame.h"
SP 10 #include "io.h"
01d871 11 #include "gyration.h"
fedf2b 12
626811 13 ts_bool run_simulation(ts_vesicle *vesicle, ts_uint mcsweeps, ts_uint inititer, ts_uint iterations, ts_uint start_iteration){
d7a113 14     ts_uint i, j;
01d871 15     ts_double l1,l2,l3;
fe24d2 16 //     char filename[255];
d7a113 17
SP 18     centermass(vesicle);
19     cell_occupation(vesicle);
626811 20     if(start_iteration<inititer) ts_fprintf(stdout, "Starting simulation (first %d x %d MC sweeps will not be recorded on disk)\n", inititer, mcsweeps);
SP 21     for(i=start_iteration;i<inititer+iterations;i++){
d7a113 22         for(j=0;j<mcsweeps;j++){
SP 23             single_timestep(vesicle);
24         }
25         centermass(vesicle);
26         cell_occupation(vesicle);
01d871 27         gyration_eigen(vesicle, &l1, &l2, &l3);
f8e6ba 28         ts_fprintf(stdout,"Done %d out of %d iterations (x %d MC sweeps).\n",i+1,inititer+iterations,mcsweeps);
1ab449 29             dump_state(vesicle,i);
58230a 30         if(i>=inititer){
d7a113 31             write_vertex_xml_file(vesicle,i-inititer);
fe24d2 32     write_master_xml_file("test.pvd");
M 33
144784 34         //    sprintf(filename,"timestep-%05d.pov",i-inititer);
fe24d2 35         //    write_pov_file(vesicle,filename);
d7a113 36         }
SP 37     }
38     return TS_SUCCESS;
39 }
d7639a 40
SP 41 ts_bool single_timestep(ts_vesicle *vesicle){
42     ts_bool retval;
43     ts_double rnvec[3];
fedf2b 44     ts_uint i,j,b;
aec47d 45     for(i=0;i<vesicle->vlist->n;i++){
d7639a 46         rnvec[0]=drand48();
SP 47         rnvec[1]=drand48();
48         rnvec[2]=drand48();
aec47d 49         retval=single_verticle_timestep(vesicle,vesicle->vlist->vtx[i],rnvec);
d7639a 50     }
SP 51
b7c9b3 52 //    ts_int cnt=0;
fedf2b 53     for(i=0;i<3*vesicle->vlist->n;i++){
dcf17d 54 //why is rnvec needed in bondflip?
SP 55 /*        rnvec[0]=drand48();
56         rnvec[1]=drand48();
57         rnvec[2]=drand48();
58 */ 
142a67 59     b=rand() % vesicle->blist->n;
d7639a 60         //find a bond and return a pointer to a bond...
SP 61         //call single_bondflip_timestep...
142a67 62         retval=single_bondflip_timestep(vesicle,vesicle->blist->bond[b],rnvec);
b7c9b3 63 //    if(retval==TS_SUCCESS) cnt++;        
fedf2b 64     }
M 65
66     for(i=0;i<vesicle->poly_list->n;i++){
58230a 67         for(j=0;j<vesicle->poly_list->poly[i]->vlist->n;j++){
M 68             rnvec[0]=drand48();
69             rnvec[1]=drand48();
70             rnvec[2]=drand48();
71             retval=single_poly_vertex_move(vesicle,vesicle->poly_list->poly[i],vesicle->poly_list->poly[i]->vlist->vtx[j],rnvec);    
72         }
fedf2b 73     }
M 74
58230a 75
M 76     for(i=0;i<vesicle->filament_list->n;i++){
77         for(j=0;j<vesicle->filament_list->poly[i]->vlist->n;j++){
78             rnvec[0]=drand48();
79             rnvec[1]=drand48();
80             rnvec[2]=drand48();
81             retval=single_filament_vertex_move(vesicle,vesicle->filament_list->poly[i],vesicle->filament_list->poly[i]->vlist->vtx[j],rnvec);    
82         }
fedf2b 83     }
M 84  
58230a 85
fedf2b 86 //    printf("Bondflip success rate in one sweep: %d/%d=%e\n", cnt,3*vesicle->blist->n,(double)cnt/(double)vesicle->blist->n/3.0);
41a035 87     if(retval);
d7639a 88     return TS_SUCCESS;
SP 89 }
90
91
92