I am trying to write a general program that involves a general (or user-defined) function h() which outputs a vector-valued function or its derivative (which is a matrix). I need to run this function on every row of a matrix X. Currently, using apply(), I get an output matrix where the outputs of h() are columns of the matrix (the derivative of h(), which is a matrix, are flattened out into vectors). I like the structure of this output of apply(), it is good for my next calculation. However, it is too slow, I wonder if there are faster ways to do this?

I have an n by m matrix X, and an array A of dimension (S,m,p). I used want to get an array B of dimension (S,n,p) where B[i,,] equals to X %*% A[i,,].