It seems the hyperparameters / architecture is not working for this kind of data.
I could overfit a small example using a lower learning rate:
data = torch.randn(1, 1, 5, 10, 10)
label = torch.zeros_like(data)
label[0, 0, 2:4, 2:7, 2:7] = 1
label.unique()
data = data.to(device)
label = label.to(device)

model = nn.Sequential(
        nn.Conv3d(1, 32, 3, 1, 1),
        nn.ReLU(),
        nn.Conv3d(32, 1, 3, 1, 1),
        nn.ReLU(),
        nn.Conv3d(1, 1, 3, 1, 1, bias=False)
    )

model = model.to(device)

loss_criterion = nn.BCEWithLogitsLoss()

optimizer = torch.optim.Adam(model.parameters(), lr=1e-5)

for epoch in range(800):
    optimizer.zero_grad()

    output = model(data)
    loss = loss_criterion(output, label)
    loss.backward()
    optimizer.step()

    print('Epoch {}, loss {}'.format(epoch, loss.item()))

print('training finished')


pred = torch.sigmoid(output) > 0.5
print((pred.float() == label).float().sum() / label.nelement())

If you run your code, check the bias of your last conv layer, as in my runs apparently the activation was “killed” and the model output was a constant value given by the bias.
In my experience this happens if your learning rate is too high and throws the parameters off.
PS: The threshold of 0.5 was of course a bad idea, as I didn’t realize you are using logits, not probabilities.