Create a map sized CxWxH of the feature a_i repeated WxH times, then calculate a dot product along the channel dimension between that and your original feature map of size CxWxH (let’s call it A):
A_rep = a_i.repeat((1,W,H))
result = (A_rep*A).sum(axis=0)

I ignored N for simplicity but the concept is similar with N>1