Calculate/approximate the velocity and rotation impulse on a rigid body, given only a force at a position and the rigidbody's center of mass?

Question

I'm working on a compute shader that allows two-way physics interaction between a fluid particle engine (NVidia Flex) and a rigidbody physics engine (Unity3D Engine).

Basically for the compute shader I'm using 5 input buffers:
(float4) particle velocities / shape contact indices
(float4) shape centroids
(int) shape flags - dynamic vs static etc
(int) particle indices
(float4) particle positions and mass

and 2 output buffers:
(float3) velocity delta
(float3) rotational velocity delta

The functionality that I'm looking for is minimal and does not need to be accurate, it just needs to be somewhat believable, as I'm mostly just using this for visual effects. I know I can create rigidbody constraints with NVidia flex particles and use that, but in my case this would not be practical because my fluid simulation uses very small particles and medium-sized rigiddodies would use many more particles with the NVidia rigid constraints than the documentation says is recommended per body.

So anyway, I've gotten to the point where in my shader all I need is a physics formula to take in the origin point of a force in world space, the force vector, the shape's center of mass in world space, and I need it to give me both the net delta velocity of the shape (assuming uniform density), and the net rotational velocity of the shape. This function will be applied on each shape many times for each contact between itself and a particle.

Here is some psuedo code:

// The velocity of the particle at the time of contact
float4 contactVelocity;

// The index of the shape that the particle collided with
int shapeIndex;

// The particle's position in the world (which should be the same as the contact point)
float3 pos;

// The mass of the particle
float mass;

// The shape's center of mass
float3 shapeOrigin;

// TODO: define ApplyVelForce & ApplyRotForce
velDelta[shapeIndex] = velDelta[shapeIndex] + GetVelForce(shapeOrigin, contactVelocity * mass, pos);
rotVelDelta[shapeIndex] = rotVelDelta[shapeIndex] + GetRotForce(shapeOrigin, contactVelocity * mass, pos);

// function definitions
float3 GetVelForce(float3 shapeCentroid, float3 force, float3 forcePoint){ /* TODO define */ }
float3 GetRotForce(float3 shapeCentroid, float3 force, float3 forcePoint){ /* TODO define */ }

If anyone knows a relatively simple formula to calculate or even approximate these velocity and rotation forces reasonably efficiently, please let me know. I've scoured google but all the articles about this seem to be way over my head. I just don't think I've got enough experience and knowledge about kinematics yet to figure out formulas like these on my own.

Answer 1

What you want is in general underdetermined, given the information you want as input. In the general case, you need as input:

Input:

center of mass of shape : x_c (changes with time)
orientation of shape in space : U (orthogonal matrix, equivalent to a coordinate frame
                                    attached to and moving with the shape, 
                                    changes with time)
Force vector in world space : f (changes with time)
point of force in world space : x_F (changes with time)
inertia matrix of the shape : J (constant) 
mass of shape : m (constant)

Output:
Velocity of center of mass : v_c (changes with time)
Angular velocity of shape in world frame : o (changes with time)

In general, the full description of a rigid body's position and orientation in world frame at any given moment of time is given by:

x_c = x_c(t), U = U(t)

where U is 3 by 3 orthogonal matrix, whose columns are the coordinates of 3 unit vectors of a coordinate frame firmly attached to and moving with the rigid body with origin at the point x_c.

Assume that some force is applied to a point X_F on the rigid body, where the coordinates of X_F are in the frame attached to the rigid body (not in the world system). In reasonable models we assume that he point X_F doesn't change with time in the body-fixed frame (of course it changes in the world frame). Let the vector of the force, in the frame attached to the body, be F and the angular velocity expressed in the body-fixed frame be O.

Then, in world frame:

point of force: 
x_F = x_c + U * X_F

force:
f = U*F

angular velocity (along the momentary axis of rotation):
o = U*O

The equations of motion then are

d/dt(x_c) = v_c

d/dt(v_c) = U*F / m

d/dt (J*O) = cross(J*O, O) + U.T * cross(X_F, F)

d/dt U = U*matrix(O)

where

x_c = center of mass,  
v_c = velocity of center of mass,  
O = angular velocity in body fixed frame,  
U = 3 by 3 orthogonal a matrix with columns the 
    orientations of the axis of the body fixed system 
F = force in body-fixed frame
X_F = point of force in body fixed frame

So, a simple iterative scheme that produces the motion of a rigid body is

input:
x_c_in, v_c_in, U_in, O_in, J, m, X_F, F, h = time-step 

x_c_out = x_c_in + h*v_c_in
v_c_out = v_c_in + h*(U_in * F_in)/m
Rot_hO = rotate(h*O_in) 
        (perform rotation around axis determined by vector O and angle h*norm(O))
U_out = U_in * Rot_hO
O_out = O_in + h*inverse(J)*( cross(J*O_in, O_in) + cross(X_F, F) )

output:
x_c_out, v_c_out, U_out, O_out,

Somewhere there should be a function that calculates the force F in the body-fixed frame. It could be that F is simply a constant vector, it could be changing its magnitude and direction only with respect to time F = F(t), or it could in general depend on the body position, orientation, velocity and even angular velocity, so F = F(x_c, v_c, U, O).