 __device__ inline real4 DifVelocityRho(
 		const real3 & dist3,
 		const real_ & d,
-		const real_ & rSPH,
 		const real4 & velMasA,
 		const real3 & vel_XSPH_A,
 		const real4 & velMasB,

//					}
					real4 derivVelRho = R4(0.0f);
					real3 vel_XSPH_B = FETCH(vel_XSPH_Sorted_D, j);
					derivVelRho = DifVelocityRho(dist3, d, rSPH, velMasA, vel_XSPH_A, velMasB, vel_XSPH_B, rhoPresMuA, rhoPresMuB, multViscosit);
					derivV += R3(derivVelRho);
					derivRho += derivVelRho.w;
				}
				else if (fabs(rhoPresMuA.w - rhoPresMuB.w) > 0) { //implies: one of them is solid/boundary, ther other one is solid/boundary of different type or different solid
//					real3 dV = DifVelocity_SSI_DEM(dist3, d, rSPH, velMasA, velMasB);
					real3 dV = DifVelocity_SSI_Lubrication(dist3, d, rSPH, velMasA, velMasB);


					if (rhoPresMuA.w > 0 && rhoPresMuA.w <= numObjectsD.numRigidBodies) { //i.e. rigid

__device__ inline real3 DifVelocity_SSI_DEM(
				const real3 & dist3,
				const real_ & d,
				const real4 & velMasA,
				const real4 & velMasB) {
//printf("** DifVelocity_SSI_DEM\n");
	real_ l = paramsD.MULT_INITSPACE * paramsD.HSML - d; //penetration distance
	if (l < 0) {
		return R3(0);
	}
