[WIP] Multinomial Logistic Regression

Continuation of the work by @larsmans
https://github.com/scikit-learn/scikit-learn/pull/2814
@larsmans I suppose this should be

    grad = np.hstack([grad, C* diff.sum(axis=0).reshape(-1, 1)])

or even better we could just do `alpha = 1./C` . and multiply it with the penalty term to make it less confusing.