+ /* virial - plus sign, as dist_ij = - dist_ji - (really??) */
+// TEST ... with a minus instead
+ ai->virial.xx-=force.x*dist_ij->x;
+ ai->virial.yy-=force.y*dist_ij->y;
+ ai->virial.zz-=force.z*dist_ij->z;
+ ai->virial.xy-=force.x*dist_ij->y;
+ ai->virial.xz-=force.x*dist_ij->z;
+ ai->virial.yz-=force.y*dist_ij->z;
+
+#ifdef DEBUG
+if(ai==&(moldyn->atom[0])) {
+ printf("dVji (3bp) contrib:\n");
+ printf("%f | %f\n",force.x,ai->f.x);
+ printf("%f | %f\n",force.y,ai->f.y);
+ printf("%f | %f\n",force.z,ai->f.z);
+}
+#endif
+#ifdef VDEBUG
+if(ai==&(moldyn->atom[0])) {
+ printf("dVji (3bp) contrib:\n");
+ printf("%f | %f\n",force.x*dist_ij->x,ai->virial.xx);
+ printf("%f | %f\n",force.y*dist_ij->y,ai->virial.yy);
+ printf("%f | %f\n",force.z*dist_ij->z,ai->virial.zz);
+}
+#endif
+