I ran some tests…
Notes…

Each trail is run 100 times
Each trial runs for a max of 1000 epochs because training generally stagnates before then.
Successes is the number of trial runs that succeeded in solving XOR
Unless otherwise stated each experiment uses the following settings

hidden_dim 5
activation between linear layers ReLU
final activation Sigmoid
loss L1Loss
learning_rate 0.03
optimizer SGD



Is MSELoss is better?
loss L1Loss
38 / 100 = 38.00% successes

loss MSELoss
78 / 100 = 78.00% successes

Conclusions: MSELoss is better for this model and this data.
Is Adam(lr=0.001) is better than SGD(lr=0.03)?
loss L1Loss
learning_rate 0.001
optimizer Adam
22 / 100 = 22.00% successes
2 / 100 = 2.00% all zero predictions

loss MSELoss
learning_rate 0.001
optimizer Adam
16 / 100 = 16.00% successes

Conclusion: For L1Loss and for MSELoss Adam(lr=0.001) is quite a lot worse than SGD(lr=0.03).
What about using Sigmoid instead of ReLU?
Motivation: ReLU is piecewise linear, so maybe it doesn’t help that much.
activation Sigmoid
loss L1Loss
18 / 100 = 18.00% successes

activation Sigmoid
loss MSELoss
87 / 100 = 87.00% successes

Conclusion: Sigmoid activation between the linear layers improves the model when using MSELoss, but worsens the model when using L1Loss.
What about ELU?
activation ELU
loss L1Loss
18 / 100 = 18.00% successes

activation ELU
loss MSELoss
100 / 100 = 100.00% successes

Wow! ELU activation between the layers is really beneficial when using MSELoss, but somewhat harmful when using L1Loss.
Any other ideas?
For reference: here is my ugly/hacky code
import numpy as np

import argparse

import torch
from torch import nn
from torch.autograd import Variable
from torch import FloatTensor
from torch import optim

use_cuda = torch.cuda.is_available()

X = xor_input = np.array([[0,0], [0,1], [1,0], [1,1]])
Y = xor_output = np.array([[0,1,1,0]]).T

# Converting the X to PyTorch-able data structure.
X_pt = Variable(FloatTensor(X))
X_pt = X_pt.cuda() if use_cuda else X_pt

# Converting the Y to PyTorch-able data structure.
Y_pt = Variable(FloatTensor(Y), requires_grad=False)
Y_pt = Y_pt.cuda() if use_cuda else Y_pt

input_dim = 2
output_dim = 1

parser = argparse.ArgumentParser(description='PyTorch example')
parser.add_argument('--loss', type=str, default="L1Loss")
parser.add_argument('--hidden_dim', type=int, default=5)
parser.add_argument('--learning_rate', type=float, default=0.03)
parser.add_argument('--activation', type=str, default="ReLU")
parser.add_argument('--optimizer', type=str, default="SGD")
args = parser.parse_args()
print(args)

hidden_dim = args.hidden_dim
criterion = getattr(nn, args.loss)()
learning_rate = 0.03
activation_between_linear_layers = getattr(nn, args.activation)
optimizer_class = getattr(optim, args.optimizer)

num_epochs = 1000
num_trials = 100
count_successes = 0
count_all_zeros = 0

try:
    for trial in range(num_trials):
        print("Trial", trial + 1, end="\r")
        model = nn.Sequential(nn.Linear(input_dim, hidden_dim),
                              activation_between_linear_layers(),
                              nn.Linear(hidden_dim, output_dim),
                              nn.Sigmoid())
        optimizer = optimizer_class(model.parameters(), lr=learning_rate)

        try:
            for epoch in range(num_epochs):
                predictions = model(X_pt)
                loss_this_epoch = criterion(predictions, Y_pt)
                if loss_this_epoch.data[0] < 1e-8:
                    break
                loss_this_epoch.backward()
                optimizer.step()
        except KeyboardInterrupt:
            raise KeyboardInterrupt()

        success = 0
        all_zero = True
        for _x, _y in zip(X_pt, Y_pt):
            prediction = model(_x)
            success += 25*(int(prediction > .5)==_y).data[0]
            all_zero = all_zero and int(prediction <= .5)
        count_successes += (success == 100)
        count_all_zeros += all_zero
except KeyboardInterrupt:
    pass
print()
print("hidden_dim", args.hidden_dim)
print("activation", args.activation)
print("loss", args.loss)
print("learning_rate", args.learning_rate)
print("optimizer", args.optimizer)
trial += 1
print("%d / %d = %.2f%% successes" % (count_successes, trial, 100.*count_successes/trial))
print("%d / %d = %.2f%% all zero predictions" % (count_all_zeros, trial, 100.*count_all_zeros/trial))