+double thermodynamic_pressure_calc(t_moldyn *moldyn) {
+
+ t_3dvec dim,*tp;
+ double u,p;
+ double scale;
+ t_atom *store;
+
+ tp=&(moldyn->tp);
+ store=malloc(moldyn->count*sizeof(t_atom));
+ if(store==NULL) {
+ printf("[moldyn] allocating store mem failed\n");
+ return -1;
+ }
+
+ /* 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);
+ 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);
+ 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);
+ 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\n",(tp->x+tp->y+tp->z)/ATM);
+
+ scale=1.0+pow(moldyn->dv/moldyn->volume,ONE_THIRD);
+
+ scale_dim(moldyn,scale,1,1,1);
+ scale_dim(moldyn,scale,1,1,1);
+ potential_force_calc(moldyn);
+
+ 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;
+
+ return sqrt(p);
+}
+