Hi,
I’m not sure, but you should check the backward method of BCE in PyTorch C source.
Here are two of them:

pytorch/Loss.cpp at master · pytorch/pytorch · GitHub 2
nn/BCECriterion.c at master · torch/nn · GitHub 2

both of them utilize d(L)/d(x) = -w (y - x) / (x - x^2) as backward method.
But they’re different in using epsilon.
First one uses d(L)/d(x) = -w*(y - x) / max((x - x^2),EPS)
and for the second one, derivative is -w*(y - x) / ((1. - x + EPS) * (x + EPS))
I go with the first one.
BTW, if you want to use GPU, you should check what the exact implementation is for it.
still I don’t think that this small difference causes your problem. It’s in the floating-point error range.