I'm afraid I've never had to work on a full-featured physics engine but I'd say it depends on how you want to manage assemblies. A first solution is to aggregate the parts into a single object updating the physical characteristics (com location, inertia tensor) indeed. However if you mention joints then I would expect the engine to simultaneously manage the parts and the additional actions from the joints, all as distinct but interacting objects.

The first method does the job for my limited needs and I can safely stick with rigid-body physics + collisions management. To deal with assemblies there is a theorem to relocate the inertia tensor: https://en.wikipedia.org/wiki/Parallel_axis_theorem (see Tensor Generalization). However for simplicity I modeled my sub parts as sums of discrete mass elements themselves, so the characteristics of any assembly can be easily obtained by summing over all the elements.

(vectors in bold, "L" subscript for linear quantities, "A" for angular, "v" for velocities, "a" for accelerations, "x" the cross product, "." the dot product)

Considering points belonging to a rigid body, the known relations are for velocities:
vL = vA x r       (eq 1)
Then if you take the cross product of both sides by r:
r x vL = r x (vA x r) = r2vA - (vA . r) r = r2vA        { for any vectors a and b, a x (b x a) = a² b - (a . b) a, and r is perpendicular to vA }
which yields:
vA = (r x vL) / r2      (eq 2)

Now for accelerations you need to differentiate { knowing d(ab)/dt = (da/dt) b + a (db/dt) }. For eq 1:
aL = aA x r + vA x vL = tangential linear acceleration + radial linear acceleration (which is not zero since the point is describing a circle)

And for eq 2:
aA = (vL x vL + r x aL) / r2 = (r x aL) / r2 

Not very nice I'm afraid! You might want to have a look at https://en.wikipedia.org/wiki/Angular_velocity

No conversions needed: the impulsed-based response model deals with both aspects and collisions modify both the linear and angular velocities of the colliding objects.

However this is only the basics because as its name says the model supposes transient events. I definitely lack the experience to properly discuss this, however it seems that a robust impulsed-based model implies an overall loop iterating over all colliding pairs and that must finally converge towards a final stable state without any contact left.

When facing the issue I switched to a "penalty-based" method, which inserts springs at the contact points (generating forces and moments to be integrated as usual). But this is no silver bullet either. A first obvious issue is how to choose and adjust the spring and damping constants. Another one is stiffness problems during the integration: the time step needs to be very tightly monitored and possibly drastically reduced to avoid sub-sampling the responses.

Anyway I'd say the best is to investigate the different methods and choose the one that best suits your needs.

If you are using wrenches for momenta or forces, then you need to use the spatial inertia matrix to convert motion twists into wrenches, 

The 6×6 spatial momentum of a particle of mass \( m \) located at \( \boldsymbol{c} \) is defined by the following block matrix

$$ \mathbf{I} = \begin{bmatrix} m & -m [\boldsymbol{c}\times] \\ m [\boldsymbol{c}\times] & - m [\boldsymbol{c}\times][\boldsymbol{c}\times] \end{bmatrix} $$

Here \( [\boldsymbol{c}\times] \) is the [3×3 matrix representation of the cross product](https://en.wikipedia.org/wiki/Cross_product#Conversion_to_matrix_multiplication).

Example:

I am using a momentum wrench because the math is easier to demonstrate, but the same rules apply to forces (since they are also wrenches). But particle accelerations are _not_ twists, You have to convert them to spatial accelerations to use them.

The velocity twist of a particle is $$ \mathbf{v} = \begin{bmatrix}\boldsymbol{v}_{C}+\boldsymbol{c}\times\boldsymbol{\omega}\\
\boldsymbol{\omega} \end{bmatrix} $$ as seen by the origin. Using the spatial inertia from above the momentum is

$$ \mathbf{p}=\mathbf{I}\mathbf{v} $$ which expands out to

$$ \begin{bmatrix}\boldsymbol{p}\\
\boldsymbol{c}\times\boldsymbol{p}
\end{bmatrix}=\begin{bmatrix}m & -m[\boldsymbol{c}\times]\\
m[\boldsymbol{c}\times] & -m[\boldsymbol{c}\times][\boldsymbol{c}\times]
\end{bmatrix}\begin{bmatrix}\boldsymbol{v}_{C}+\boldsymbol{c}\times\boldsymbol{\omega}\\
\boldsymbol{\omega}
\end{bmatrix} $$

carry out the multiplication

$$\begin{bmatrix}m\boldsymbol{v}_{C}+m\left(\boldsymbol{c}\times\boldsymbol{\omega}\right)-m\left(\boldsymbol{c}\times\boldsymbol{\omega}\right)\\
m\boldsymbol{c}\times\left(\boldsymbol{v}_{C}+\boldsymbol{c}\times\boldsymbol{\omega}\right)-m\boldsymbol{c}\times\left(\boldsymbol{c}\times\boldsymbol{\omega}\right)
\end{bmatrix} = \begin{bmatrix}m\boldsymbol{v}_{C}\\
\boldsymbol{c}\times\left(m\boldsymbol{v}_{C}\right)
\end{bmatrix}$$

The above confirms the known vector relationships of \( \boldsymbol{p} = m \boldsymbol{v}_C \) and \( \boldsymbol{L} = \boldsymbol{c} \times \boldsymbol{p} \)