MS_INTERPOLATE - Symmetry preserving elasticity interpolation

Symmetry aware linear interpolation between two elasticity matrices.

[ Cint, RhoInt ] = MS_interpolate( C1, rho1, C2, rho2, frac)

Usage:
    [ Cint, RhoInt ] = MS_interpolate( C1, rho1, C2, rho2, frac)
        Returns an elasticity matrix interpolated frac distance
        between elasticity matricies C1 and C2. Frac = 1 gives
        the elasticity of C1, frac = 0 gives the elasticity of
        C2.

Notes:
    A common approach to the interpolation of elasticity matricies
    is to interpolate elements component-wise. This fails for cases
    where the elasticity matricies are just rotated versions of each
    other. The approach taken here is to first rotate the two matricies
    into a common orientation (the one generated by MS_axes is chosen)
    before interpolating element-wise (Voigt-Reuss-Hill is used). The
    two rotations into the common orientation are then interpolated
    and the interpolated matrix rotated back the the starting
    orientation. Shoemake's (1985) slerp approach is used to interpolate
    the orientations. Rotation matricies are converted to quaternions
    (and back again) using John Fuller's SpinCalc function.

References:
    Shoemake, K (1985) "Animating rotation with quaternion curves" ACM
    SIGGRAPH Computer Graphics vol.19 pp.245-254.

See also: MS_AXES, MS_VRH