I can’t recommend any specific library for hyperparameter optimization, as I’m usually not using them, so let’s wait for others to chime in.
Using your model and removing the last activation I’m able to overfit a random input and target tensor:
class Net(nn.Module):
    def __init__(self, n_features):
        super(Net, self).__init__()
        # n_features = no. inputs
        n1 = 8 # no. nodes in layer 1
        n2 = 5 # no. nodes in layer 2
        n3 = 4 # no. nodes in layer 3
        n4 = 5 # no. nodes in layer 4
        n5 = 2 # no. nodes in layer 5
        self.fc1 = nn.Linear(n_features,n1)
        self.fc2 = nn.Linear(n1,n2)
        self.fc3 = nn.Linear(n2,n3)
        self.fc4 = nn.Linear(n3,n4)
        self.fc5 = nn.Linear(n4,n5)
        self.fc6 = nn.Linear(n5,1)
    def forward(self, x):
        x = torch.tanh(self.fc1(x))
        x = torch.sigmoid(self.fc2(x)) 
        x = torch.sigmoid(self.fc3(x)) 
        x = torch.tanh(self.fc4(x)) 
        x = torch.tanh(self.fc5(x))
        return self.fc6(x)

torch.manual_seed(2809)
X_train = torch.randn(64, 10)
y_train = torch.randn(64, 1)
net = Net(X_train.shape[1])
criterion = nn.MSELoss() # loss function
optimizer = torch.optim.Adam(net.parameters(), lr=1e-3)

## training 
for epoch in range(15000):
    y_pred = net(X_train)
    train_loss = criterion(y_pred, y_train)
    optimizer.zero_grad()        
    train_loss.backward()
    optimizer.step()
    print('epoch {}, loss {}'.format(epoch, train_loss.item()))

plt.figure()
plt.plot(y_pred.detach().numpy(), label='predicted')
plt.plot(y_train.numpy(), label='actual')
plt.ylabel('output y')
plt.legend()
plt.show()


You should also check, why you are using squeeze() on the model output, as I could imagine it might create unwanted broadcasting in the loss calculation (the criterion should raise a warning).