Single index contraction in sparse matrix multiplication of multidimensional array

I’d like to be able to effectively do something like C[a,b,c] = A[b,d]*B[a,d,c], where A is a sparse matrix, and C and B are dense (and C and B could be other ranks, including 1). This is indeed just contraction of a single index, but the catch is that the contraction index is not know until runtime. I briefly looked at some tensor/Einstein sum packages, but can’t seem to find anything that would be able to do this in-place (possible I didn’t look hard enough!).

I’m currently thinking of doing something along the lines of @views C[indices...] = A * B[indices...], where indices is a tuple containing a single Colon(). I’m not sure if doing this with @views and splatting destroys performance due to non-contiguous memory layout. Either way, I can’t figure out if there is an efficient way to construct the indices tuple by iterating over all of the other (non-contracted) indices when the rank of B/C is not known until runtime. Is this possible without effectively writing different routines for different ranks? I could always change my data structure to not use multidimensional arrays, but that would likely complicate the later sparse matrix construction, and I’d like to avoid that if possible.

It’s quite likely that I’m doing something stupid, and missing the point entirely, but if there is a way to do this concisely and without a crippling hit to performance it would help me a lot with the prototyping stage I’m currently in with Julia.