Thanks! I think that did the trick, i.e. registering a hook function for each matrix:
def _E_hook(grad):
    return grad * Variable(self.E_mask)

def _D_hook(grad):
    return grad * Variable(self.D_mask)

Just need to convince myself that computing the gradients of, e.g.
output w.r.t D and then zeroing out the D gradients before backpropogating is equivalent to never computing gradients for certain dimensions of D in the first place.