diff --git a/src/ComputeGridGrid.cu b/src/ComputeGridGrid.cu index f6cb72cc78899b85ca0395e28e9021dedd9c92ed..f2eb6fde0bd9802b89537f39af10efd9f32618df 100644 --- a/src/ComputeGridGrid.cu +++ b/src/ComputeGridGrid.cu @@ -2,7 +2,6 @@ #include "ComputeGridGrid.cuh" #include "RigidBodyGrid.h" #include "CudaUtil.cuh" -//RBTODO handle periodic boundaries //RBTODO: add __restrict__, benchmark (Q: how to restrict member data?) class GridPositionTransformer { @@ -17,6 +16,7 @@ private: const Vector3 c; const BaseGrid* s; }; + //class PmfPositionTransformer : public BasePositionTransformer { class PmfPositionTransformer { public: @@ -52,8 +52,7 @@ inline void common_computeGridGridForce(const RigidBodyGrid* rho, const RigidBod Vector3 r_pos= rho->getPosition(r_id); /* i,j,k value of voxel */ r_pos = basis_rho.transform( r_pos ); - r_pos = transformer(r_pos); - const Vector3 u_ijk_float = basis_u_inv.transform( r_pos ); + const Vector3 u_ijk_float = basis_u_inv.transform( transformer( r_pos ) ); // RBTODO: Test for non-unit delta /* Vector3 tmpf = Vector3(0.0f); */ /* float tmpe = 0.0f; */ @@ -71,7 +70,8 @@ inline void common_computeGridGridForce(const RigidBodyGrid* rho, const RigidBod const float r_val = rho->val[r_id]; /* maybe move to beginning of function? */ force[tid].f = basis_u_inv.transpose().transform( r_val*(force[tid].f) ); /* transform to lab frame, with correct scaling factor */ force[tid].e = r_val; - // Calculate torque about origin_u in the lab frame + + // Calculate torque about origin_rho in the lab frame torque[tid].f = r_pos.cross(force[tid].f); } diff --git a/src/RigidBodyController.cu b/src/RigidBodyController.cu index b7fd3c8964273b336705e4f912e1f7d4df8eedd3..13438b8c758384996f3c72309d2bc62071a1f7c7 100644 --- a/src/RigidBodyController.cu +++ b/src/RigidBodyController.cu @@ -902,10 +902,10 @@ void RigidBodyForcePair::processGPUForces(BaseGrid* sys) { tmpT = tmpT + torques[i][j]; } - // tmpT is the torque calculated about the origin of grid k2 (e.g. c2) + // tmpT is the torque calculated about the origin of density grid // so here we transform torque to be about rb1 - Vector3 o2 = getOrigin2(i); - tmpT = tmpT - sys->wrapDiff(rb1->getPosition() - o2).cross( tmpF.f ); + Vector3 o1 = getOrigin1(i); + tmpT = tmpT - (rb1->getPosition() - o1).cross( tmpF.f ); // clear forces on GPU gpuErrchk(cudaMemset((void*)(forces_d[i]),0,nb*sizeof(ForceEnergy)));