Constant volume works. CAVEAT: in centermass I had to recalculate all normals of triangles again! This could influence behaviour of the system.
| | |
| | | } |
| | | // All checks OK! |
| | | |
| | | for(j=0;j<vtx_moved->neigh_no;j++){ |
| | | memcpy((void *)&backupvtx[j+1],(void *)vtx_moved->neigh[j],sizeof(ts_vertex)); |
| | | } |
| | | dvol=0.0; |
| | | for(j=0;j<vtx_moved->tristar_no;j++){ |
| | | dvol-=vtx_moved->tristar[j]->volume; |
| | |
| | | vtx_moved->z=vtx_moved->z*(1-dh/Rv); |
| | | //check for constraints |
| | | if(constvolConstraintCheck(vesicle, vtx_moved)==TS_FAIL){ |
| | | for(j=0;j<vtx_moved->neigh_no;j++){ |
| | | memcpy((void *)vtx_moved->neigh[j],(void *)&backupvtx[j+1],sizeof(ts_vertex)); |
| | | } |
| | | vtx_moved=memcpy((void *)vtx_moved,(void *)&backupvtx[0],sizeof(ts_vertex)); |
| | | //also, restore normals |
| | | for(j=0;j<vtx_moved->tristar_no;j++) triangle_normal_vector(vtx_moved->tristar[j]); |
| | |
| | | if(fabs(voldiff)/vesicle->volume < vesicle->tape->constvolprecision){ |
| | | //calculate energy, return change in energy... |
| | | // fprintf(stderr, "Constvol success! %e\n",voldiff); |
| | | for(j=0;j<vtx_moved->neigh_no;j++){ |
| | | memcpy((void *)&backupvtx[j+1],(void *)vtx_moved->neigh[j],sizeof(ts_vertex)); |
| | | } |
| | | |
| | | oenergy=vtx_moved->energy; |
| | | energy_vertex(vtx_moved); |
| | | delta_energy=vtx_moved->xk*(vtx_moved->energy - oenergy); |
| | |
| | | ts_bool constvolumerestore(ts_vertex *vtx_moved,ts_vertex *vtx_backup){ |
| | | ts_uint j; |
| | | memcpy((void *)vtx_moved,(void *)&vtx_backup[0],sizeof(ts_vertex)); |
| | | for(j=0;j<vtx_moved->neigh_no;j++){ |
| | | memcpy((void *)vtx_moved->neigh[j],(void *)&vtx_backup[j+1],sizeof(ts_vertex)); |
| | | } |
| | | for(j=0;j<vtx_moved->tristar_no;j++) triangle_normal_vector(vtx_moved->tristar[j]); |
| | | for(j=0;j<vtx_moved->neigh_no;j++){ |
| | | // memcpy((void *)vtx_moved->neigh[j],(void *)&vtx_backup[j+1],sizeof(ts_vertex)); |
| | | energy_vertex(vtx_moved->neigh[j]); |
| | | } |
| | | |
| | | free(vtx_backup); |
| | | return TS_SUCCESS; |
| | |
| | | #include "cell.h" |
| | | #include "frame.h" |
| | | |
| | | |
| | | #include "triangle.h" |
| | | ts_bool centermass(ts_vesicle *vesicle){ |
| | | ts_uint i,j, n=vesicle->vlist->n; |
| | | ts_vertex **vtx=vesicle->vlist->vtx; |
| | |
| | | vesicle->cm[1]=0.0; |
| | | vesicle->cm[2]=0.0; |
| | | |
| | | for(i=0;i<vesicle->tlist->n;i++){ |
| | | triangle_normal_vector(vesicle->tlist->tria[i]); |
| | | } |
| | | |
| | | |
| | | return TS_SUCCESS; |
| | | } |
| | | |
| | |
| | | pressure=0.0 |
| | | |
| | | #Constant volume constraint (0 disable constant volume, 1 enable) |
| | | constvolswitch=1 |
| | | constvolswitch=0 |
| | | constvolprecision=1e-14 |
| | | |
| | | |
| | |
| | | for(i=start_iteration;i<inititer+iterations;i++){ |
| | | vmsr=0.0; |
| | | bfsr=0.0; |
| | | /* vesicle_volume(vesicle); |
| | | fprintf(stderr,"Volume before TS=%1.16e\n", vesicle->volume); */ |
| | | for(j=0;j<mcsweeps;j++){ |
| | | single_timestep(vesicle, &vmsrt, &bfsrt); |
| | | vmsr+=vmsrt; |
| | | bfsr+=bfsrt; |
| | | } |
| | | /* |
| | | vesicle_volume(vesicle); |
| | | fprintf(stderr,"Volume after TS=%1.16e\n", vesicle->volume); */ |
| | | vmsr/=(ts_double)mcsweeps; |
| | | bfsr/=(ts_double)mcsweeps; |
| | | centermass(vesicle); |
| | |
| | | } |
| | | |
| | | ts_bool single_timestep(ts_vesicle *vesicle,ts_double *vmsr, ts_double *bfsr){ |
| | | // vesicle_volume(vesicle); |
| | | // fprintf(stderr,"Volume before TS=%1.16e\n", vesicle->volume); |
| | | ts_bool retval; |
| | | ts_double rnvec[3]; |
| | | ts_uint i,j, b; |
| | |
| | | //find a bond and return a pointer to a bond... |
| | | //call single_bondflip_timestep... |
| | | retval=single_bondflip_timestep(vesicle,vesicle->blist->bond[b],rnvec); |
| | | // b++; retval=TS_FAIL; |
| | | if(retval==TS_SUCCESS) bfsrcnt++; |
| | | } |
| | | |
| | |
| | | // printf("Bondflip success rate in one sweep: %d/%d=%e\n", cnt,3*vesicle->blist->n,(double)cnt/(double)vesicle->blist->n/3.0); |
| | | *vmsr=(ts_double)vmsrcnt/(ts_double)vesicle->vlist->n; |
| | | *bfsr=(ts_double)bfsrcnt/(ts_double)vesicle->vlist->n/3.0; |
| | | // vesicle_volume(vesicle); |
| | | // fprintf(stderr,"Volume after TS=%1.16e\n", vesicle->volume); |
| | | return TS_SUCCESS; |
| | | } |
| | | |