function c1 = c1_calc_alt(H_plus,p)

K_a1 = p(9);
K_a2 = p(10);
K_w =  p(11);

c1 = (K_a1/H_plus + (2*K_a1*K_a2)/H_plus^2)/...
    (1+K_a1/H_plus+K_a1*K_a2/H_plus^2);

end