+ /* save unscaled potential energy + atom/dim configuration */
+ u=moldyn->energy;
+ memcpy(store,moldyn->atom,moldyn->count*sizeof(t_atom));
+ dim=moldyn->dim;
+
+ /* derivative with respect to x direction */
+ scale=1.0+moldyn->dv/(moldyn->dim.y*moldyn->dim.z);
+ scale_dim(moldyn,scale,TRUE,0,0);
+ scale_atoms(moldyn,scale,TRUE,0,0);
+ link_cell_shutdown(moldyn);
+ link_cell_init(moldyn,QUIET);
+ potential_force_calc(moldyn);
+ tp->x=(moldyn->energy-u)/moldyn->dv;
+ p=tp->x*tp->x;
+
+ /* restore atomic configuration + dim */
+ memcpy(moldyn->atom,store,moldyn->count*sizeof(t_atom));
+ moldyn->dim=dim;
+
+ /* derivative with respect to y direction */
+ scale=1.0+moldyn->dv/(moldyn->dim.x*moldyn->dim.z);
+ scale_dim(moldyn,scale,0,TRUE,0);
+ scale_atoms(moldyn,scale,0,TRUE,0);
+ link_cell_shutdown(moldyn);
+ link_cell_init(moldyn,QUIET);
+ potential_force_calc(moldyn);
+ tp->y=(moldyn->energy-u)/moldyn->dv;
+ p+=tp->y*tp->y;
+
+ /* restore atomic configuration + dim */
+ memcpy(moldyn->atom,store,moldyn->count*sizeof(t_atom));
+ moldyn->dim=dim;
+
+ /* derivative with respect to z direction */
+ scale=1.0+moldyn->dv/(moldyn->dim.x*moldyn->dim.y);
+ scale_dim(moldyn,scale,0,0,TRUE);
+ scale_atoms(moldyn,scale,0,0,TRUE);
+ link_cell_shutdown(moldyn);
+ link_cell_init(moldyn,QUIET);
+ potential_force_calc(moldyn);
+ tp->z=(moldyn->energy-u)/moldyn->dv;
+ p+=tp->z*tp->z;
+
+ /* restore atomic configuration + dim */
+ memcpy(moldyn->atom,store,moldyn->count*sizeof(t_atom));
+ moldyn->dim=dim;
+
+ printf("dU/dV komp addiert = %f %f %f\n",tp->x,tp->y,tp->z);
+
+ scale=1.0+pow(moldyn->dv/moldyn->volume,ONE_THIRD);
+
+printf("debug: %f %f\n",moldyn->atom[0].r.x,moldyn->dim.x);
+ scale_dim(moldyn,scale,1,1,1);
+ scale_atoms(moldyn,scale,1,1,1);
+ link_cell_shutdown(moldyn);
+ link_cell_init(moldyn,QUIET);
+ potential_force_calc(moldyn);
+printf("debug: %f %f\n",moldyn->atom[0].r.x,moldyn->dim.x);
+
+ printf("dU/dV einfach = %f\n",((moldyn->energy-u)/moldyn->dv)/ATM);
+
+ /* restore atomic configuration + dim */
+ memcpy(moldyn->atom,store,moldyn->count*sizeof(t_atom));
+ moldyn->dim=dim;
+
+ /* restore energy */
+ moldyn->energy=u;
+
+ link_cell_shutdown(moldyn);
+ link_cell_init(moldyn,QUIET);
+
+ return sqrt(p);
+}
+
+double get_pressure(t_moldyn *moldyn) {
+
+ return moldyn->p;
+
+}
+
+int scale_dim(t_moldyn *moldyn,double scale,u8 x,u8 y,u8 z) {
+
+ t_3dvec *dim;
+
+ dim=&(moldyn->dim);
+
+ if(x) dim->x*=scale;
+ if(y) dim->y*=scale;
+ if(z) dim->z*=scale;