Trisurf Monte Carlo simulator
Samo Penic
2012-02-23 e016c401c98b6a33392b551b2c40fccded21efe3
src/energy.c
@@ -13,6 +13,7 @@
    for(i=0;i<vlist->n;i++){
        energy_vertex(vtx[i]);
    }
    return TS_SUCCESS;
@@ -38,6 +39,7 @@
        j=data->neigh[jj-1];
        jp=data->neigh[jjp-1];
        jm=data->neigh[jjm-1];
//        printf("tristar_no=%u, neigh_no=%u, jj=%u\n",data->tristar_no,data->neigh_no,jj);
        jt=data->tristar[jj-1];
        x1=vtx_distance_sq(vtx,jp); //shouldn't be zero!
        x2=vtx_distance_sq(j,jp); // shouldn't be zero!
@@ -72,24 +74,24 @@
        tot=0.5*tot;
        xlen=vtx_distance_sq(j,vtx);
#ifdef  TS_DOUBLE_DOUBLE 
        data->bond[jj-1]->data->bond_length=sqrt(xlen);
        data->bond[jj-1]->bond_length=sqrt(xlen);
#endif
#ifdef  TS_DOUBLE_FLOAT
        data->bond[jj-1]->data->bond_length=sqrtf(xlen);
        data->bond[jj-1]->bond_length=sqrtf(xlen);
#endif
#ifdef  TS_DOUBLE_LONGDOUBLE 
        data->bond[jj-1]->data->bond_length=sqrtl(xlen);
        data->bond[jj-1]->bond_length=sqrtl(xlen);
#endif
        data->bond[jj-1]->data->bond_length_dual=tot*data->bond[jj-1]->data->bond_length;
        data->bond[jj-1]->bond_length_dual=tot*data->bond[jj-1]->bond_length;
        s+=tot*xlen;
        xh+=tot*(j->data->x - data->x);
        yh+=tot*(j->data->y - data->y);
        zh+=tot*(j->data->z - data->z);
        txn+=jt->data->xnorm;
        tyn+=jt->data->ynorm;
        tzn+=jt->data->znorm;
        txn+=jt->xnorm;
        tyn+=jt->ynorm;
        tzn+=jt->znorm;
    }
    
    h=xh*xh+yh*yh+zh*zh;
@@ -116,7 +118,9 @@
        data->curvature=-sqrtl(h);
    }
#endif
//TODO: MAJOR!!!! What is vtx->data->c?????????????? Here it is 0!
// What is vtx->data->c?????????????? Here it is 0!
// c is forced curvature energy for each vertex. Should be set to zero for
// norman circumstances.
    data->energy=0.5*s*(data->curvature/s-data->c)*(data->curvature/s-data->c);
    return TS_SUCCESS;