Hey,
Here is an implementation that should be significantly faster and match what you want. Because of the very different low precision of your ground truth, the test trips, but comparing the two forward passes lead to very small error (same up to 5 places):
class Conv2D(torch.nn.Module):
    def __init__(self, in_channels, out_channels, kernel_size, stride=1, padding=0):
        super(Conv2D, self).__init__()
        # These are ones only for testing
        self.weight = Parameter(torch.ones(out_channels, in_channels, kernel_size, kernel_size), requires_grad=False)
        self.bias = Parameter(torch.ones(out_channels), requires_grad=False)
        self.kernel_size = kernel_size
        self.stride = stride
        self.padding = padding

    def forward(self, img):
        # img is (num_channels x n x n x num_errors)
        # outer dimension of each image channel (assumes n x n)
        num_channels = img.size(0)
        n = img.shape[-2]
        # kernel is square
        k = self.kernel_size
        # padding is the same in both H and W
        p = self.padding

        # Consider the last dimension as a batch
        img = img.permute(-1, 0, 1, 2)
        # Use regular conv
        out = F.conv2d(img, self.weight, None, self.stride, self.padding)
        # Add custom biases only on the first entry
        out[0] += self.bias.view(-1, 1, 1) * num_channels
        # Reshape to the original input
        out = out.permute(1, 2, 3, 0)
        return out

Note that you could do similar stuff with a batch of input by doing a similar permute as (0, -1, 1, 2, 3). Then collapse the first dimension in a common batch. After the conv, uncollapse the two dimensions and permute back.