Why does the choice of solver result in different weight in sklearn (LogisticRegression)

Why do different solvers result in different learned weights on simple problems in sklearn? This doesn't smell like numerical instability -- it looks like the different solvers are designed to converge on different sets of weights. What's going on?

It seems to be related to regularization, as the effect goes away with increasing C

%matplotlib inline
import time
import pandas as pd
import numpy as np
from sklearn.linear_model.logistic import LogisticRegression

def linear(X,b): return X.dot(b.T)
def sigmoid(X,b): return 1/(1+np.exp(-linear(X,b)))

b = np.array([1.,1.])

def mk(X,y,model):
    model.fit(X,y)
    t = time.clock()
    c = np.concatenate((model.coef_[0,:], model.intercept_[[0]]))
    return c, time.clock()-t

max_iter = 1000
tol=0.00000001
for C in [0.0001, 0.001, 0.1, 1., 10., 1000.]:
    for n in [200, 20000]:
        x = np.arange(-5., 5, 10./n)
        Xb = np.vstack([x, np.ones(len(x))]).T
        p = sigmoid(Xb,b)
        y = np.random.binomial(1, p)
        X = x.reshape((len(x), 1))
        print C, n
        for solver in ["liblinear", "newton-cg", "lbfgs", "sag"]:
            print "  {:>14}".format(solver), mk(X,y, LogisticRegression(C=C, solver=solver, penalty='l2', max_iter=max_iter, tol=tol))

The converged coefficients are below. The training data is generated with B=[1., 1.]

C      n
0.0001 200
       liblinear (array([ 0.01858615,  0.00219137]), 1.399999999995849e-05)
       newton-cg (array([ 0.01867489,  0.44809872]), 1.5000000000098268e-05)
           lbfgs (array([ 0.01867489,  0.44809872]), 1.5000000000098268e-05)
             sag (array([ 0.01867489,  0.448096  ]), 1.3000000000040757e-05)
0.0001 20000
       liblinear (array([ 0.49062599,  0.14814571]), 2.6999999999999247e-05)
       newton-cg (array([ 0.51257455,  0.59982733]), 2.6000000000081513e-05)
           lbfgs (array([ 0.51257471,  0.59982645]), 2.8999999999834714e-05)
             sag (array([ 0.51257455,  0.59982753]), 1.3000000000040757e-05)
0.001 200
       liblinear (array([ 0.15311436,  0.0154459 ]), 1.2999999999596668e-05)
       newton-cg (array([ 0.15456163,  0.34240309]), 1.3000000000040757e-05)
           lbfgs (array([ 0.15456166,  0.34240284]), 1.3000000000040757e-05)
             sag (array([ 0.15456163,  0.34240308]), 1.100000000020529e-05)
0.001 20000
       liblinear (array([ 0.81980482,  0.58214727]), 3.799999999998249e-05)
       newton-cg (array([ 0.86826695,  0.87335432]), 3.0999999999892225e-05)
           lbfgs (array([ 0.86826513,  0.87335497]), 3.2000000000032e-05)
             sag (array([ 0.86826695,  0.87335433]), 2.9999999999752447e-05)
0.1 200
       liblinear (array([ 0.80161404,  0.39606419]), 2.3999999999801958e-05)
       newton-cg (array([ 0.82243805,  0.57384594]), 2.9000000000056758e-05)
           lbfgs (array([ 0.82243822,  0.57384605]), 3.0000000000196536e-05)
             sag (array([ 0.82243805,  0.57384595]), 2.4999999999941735e-05)
0.1 20000
       liblinear (array([ 0.98022283,  0.97426356]), 3.799999999998249e-05)
       newton-cg (array([ 0.98180586,  0.98060186]), 2.3999999999801958e-05)
           lbfgs (array([ 0.98181404,  0.98060174]), 3.1000000000336314e-05)
             sag (array([ 0.98180586,  0.98060185]), 3.0999999999892225e-05)
1.0 200
       liblinear (array([ 0.77546388,  0.68083769]), 2.3999999999801958e-05)
       newton-cg (array([ 0.78084612,  0.71273566]), 1.1999999999900979e-05)
           lbfgs (array([ 0.7808462 ,  0.71273586]), 1.3000000000040757e-05)
             sag (array([ 0.78084612,  0.71273566]), 1.0000000000065512e-05)
1.0 20000
       liblinear (array([ 1.01286416,  1.01699803]), 3.1000000000336314e-05)
       newton-cg (array([ 1.013046  ,  1.01769302]), 1.6000000000016e-05)
           lbfgs (array([ 1.0130471 ,  1.01769538]), 3.0999999999892225e-05)
             sag (array([ 1.013046  ,  1.01769301]), 2.0999999999826713e-05)
10.0 200
       liblinear (array([ 1.19187232,  1.52139602]), 9.999999999621423e-06)
       newton-cg (array([ 1.19737431,  1.53737059]), 4.300000000023729e-05)
           lbfgs (array([ 1.19737519,  1.53736531]), 1.8999999999991246e-05)
             sag (array([ 1.19737429,  1.53737061]), 1.700000000059987e-05)
10.0 20000
       liblinear (array([ 0.9993963,  0.9748259]), 9.19999999995369e-05)
       newton-cg (array([ 0.9994126 ,  0.97489034]), 3.0999999999892225e-05)
           lbfgs (array([ 0.99941332,  0.97489202]), 3.10000000007804e-05)
             sag (array([ 0.99941261,  0.97489035]), 1.4000000000180535e-05)
1000.0 200
       liblinear (array([ 1.31759982,  1.48775989]), 9.000000000369823e-06)
       newton-cg (array([ 1.3176617 ,  1.48792683]), 1.4000000000180535e-05)
           lbfgs (array([ 1.31766183,  1.48792601]), 4.499999999918458e-05)
             sag (array([ 1.31766164,  1.48792671]), 1.100000000064938e-05)
1000.0 20000
       liblinear (array([ 0.99608584,  0.96883584]), 4.099999999951365e-05)
       newton-cg (array([ 0.99608601,  0.96883648]), 3.300000000017178e-05)
           lbfgs (array([ 0.99608665,  0.96883797]), 2.9999999999752447e-05)
             sag (array([ 0.99608603,  0.96883652]), 1.5000000000320313e-05)

And also, ... all the versions seem to converge to different coefficients depending on population size. Its almost as if the cost is not normalized by the population size and the effective regularization is C/n. Is that the case?