I'm trying to find the value of a parameter for fit experimental data to a mode I've defined the function of the model and it works OK, but when I try to use NMinimize to perform an optimization of one of the parameters (Ca0) i obtain errors...

Pleease help...

ClearAll["Global`*"]

n = 3;                                 (*Nombre de protons*)
pKa = List[2, 7, 12];(*Llista de pKa*) 
Ka = 10^-pKa;
Kw = 1 10^-14;
(*Ca0=0.01;     *)                             (*Concentració inicial \
de l'àcid*) 
Cb0 = 0.1;                                   (*Concentració inicial \
de la base*)
Va = 50; 
Vfinal = (n + 1) Va 0.01/Cb0 ;      (*Volum inicial de àcid*)

data2 = {{0.`, 2.1951268820918415`}, {0.2`, 
   2.222186986869529`}, {0.4`, 2.2512980860016265`}, {0.6`, 
   2.2783770173176072`}, {0.8`, 2.3001099277084953`}, {1.`, 
   2.3423619242548996`}, {1.2`, 2.3547166009392333`}, {1.4`, 
   2.3996174765200364`}, {1.6`, 2.416626954675591`}, {1.8`, 
   2.4602027064867333`}, {2.`, 2.4640237846626425`}, {2.2`, 
   2.515899195269751`}, {2.4`, 2.5371201132316243`}, {2.6`, 
   2.5917325517367167`}, {2.8`, 2.635015919284972`}, {3.`, 
   2.6559037249406865`}, {3.2`, 2.7137317935286784`}, {3.4`, 
   2.768268005142373`}, {3.6`, 2.8486362809900103`}, {3.8`, 
   2.9272739626775715`}, {4.`, 2.9910766465593444`}, {4.2`, 
   3.1100059845345243`}, {4.4`, 3.2351821137310903`}, {4.6`, 
   3.3884792404928628`}, {4.8`, 3.713909504539493`}, {5.`, 
   4.657346640417365`}, {5.2`, 5.642022376185781`}, {5.4`, 
   5.920761158523168`}, {5.6`, 6.133820370113931`}, {5.8`, 
   6.277160209202424`}, {6.`, 6.396592296557952`}, {6.2`, 
   6.510348141873552`}, {6.4`, 6.58282442237076`}, {6.6`, 
   6.6841010142650354`}, {6.8`, 6.749025731626233`}, {7.`, 
   6.80537151411917`}, {7.2`, 6.882318026534458`}, {7.4`, 
   6.97411420349684`}, {7.6`, 7.033695546407676`}, {7.8`, 
   7.090825695141246`}, {8.`, 7.1697036780217065`}, {8.2`, 
   7.2541234819637594`}, {8.4`, 7.314957179725143`}, {8.6`, 
   7.419609456279908`}, {8.8`, 7.490168704985635`}, {9.`, 
   7.6185352333071545`}, {9.2`, 7.710565123853162`}, {9.4`, 
   7.875275195079028`}, {9.6`, 8.046775328312181`}, {9.8`, 
   8.365877291018675`}, {10.`, 9.311628706543383`}, {10.2`, 
   10.271309818641138`}, {10.4`, 10.574955970341009`}, {10.6`, 
   10.747384415195413`}, {10.8`, 10.869764090100094`}, {11.`, 
   10.977786342664643`}, {11.2`, 11.05945230118756`}, {11.4`, 
   11.129974570610974`}, {11.6`, 11.188505298944662`}, {11.8`, 
   11.254243597709767`}, {12.`, 11.272701078805834`}, {12.2`, 
   11.336195013103278`}, {12.4`, 11.357086120093598`}, {12.6`, 
   11.389136555905981`}, {12.8`, 11.442736581710385`}, {13.`, 
   11.468249961197614`}, {13.2`, 11.498397369523966`}, {13.4`, 
   11.534569331510793`}, {13.6`, 11.565380682448227`}, {13.8`, 
   11.573653974100603`}, {14.`, 11.60993862023272`}, {14.2`, 
   11.622773138775505`}, {14.4`, 11.661430296295993`}, {14.6`, 
   11.681742643808693`}, {14.8`, 11.678480510410138`}, {15.`, 
   11.711750459487295`}, {15.2`, 11.728743730220623`}, {15.4`, 
   11.744631262581105`}, {15.6`, 11.774701340208647`}, {15.8`, 
   11.766793922998312`}, {16.`, 11.784660551516614`}, {16.2`, 
   11.791911767914888`}, {16.4`, 11.839444506570322`}, {16.6`, 
   11.841913020633076`}, {16.8`, 11.849066998047896`}, {17.`, 
   11.868840789502134`}, {17.2`, 11.878040759973212`}, {17.4`, 
   11.89475342951758`}, {17.6`, 11.902869029005489`}, {17.8`, 
   11.921869349935383`}, {18.`, 11.925518698050961`}, {18.2`, 
   11.95565106647782`}, {18.4`, 11.935640586131818`}, {18.6`, 
   11.941835839838838`}, {18.8`, 11.952434190549537`}, {19.`, 
   11.960775787389656`}, {19.2`, 11.984355376300066`}, {19.4`, 
   11.99581142911412`}, {19.6`, 12.001987273737177`}, {19.8`, 
   12.027155634634736`}, {20.`, 12.03078325340934`}}

Vb = data2[[All, 1]];
pHexp = data2[[All, 2]];

l = Length[Vb];
pHc = {};

(*CÀLCUL DELS PUNTS DE LA VALORACIÓ*)

calcpH[Va_, Vb_List, Ca0_, Cb0_, n_, Ka_List, Kw_] := 
 Module[{Vt, H, pH}, pHc = {};
  {For[j = 1, j <= l, j++,
    {Clear[H];
     Vt = Va + Vb[[j]]; Ca = (Ca0 Va)/Vt; Cb = (Cb0 Vb[[j]])/Vt;
     sols = Solve[{(\!\(
\*UnderoverscriptBox[\(\[Sum]\), \(i = 1\), \(n\)]\((i\ \ *
\*SuperscriptBox[\(H\), \(n - i\)]\ *\(
\*UnderoverscriptBox[\(\[Product]\), \(k = 
                   1\), \(i\)]Ka[\([\)\(k\)\(]\)]\))\)\)/(H^n + \!\(
\*UnderoverscriptBox[\(\[Sum]\), \(i = 1\), \(n\)]\((
\*SuperscriptBox[\(H\), \(n - i\)]*\(
\*UnderoverscriptBox[\(\[Product]\), \(k = 
                    1\), \(i\)]Ka[\([\)\(k\)\(]\)]\))\)\))*Ca) - Cb - 
           H + Kw/H == 0 && H >= 0}, H];
     H = H /. sols // N; pH = -Log10[H];
     AppendTo[pHc, pH]}];
   Valoracio = Transpose[{Vb, Flatten[pHc]}];
   }
  ]

sseval[Ca0_] := {
  calcpH[Va, Vb, Ca0, Cb0, n, Ka, Kw];
  Print["Ca0 = ", Ca0];
  sse = Norm[pHexp - pHc]
  }

NMinimize[{sseval[Ca0], Ca0 > 0}, {Ca0}] // Quiet


titol = titol <> "\nn=" <> ToString[n] <> ", valors de pKa =" <> 
   ToString[pKa];
(*FORMAT DEL GRÀFIC*)
estil = {PlotStyle -> {Dashing[Tiny], Thickness[0.001]}, 
   ImageSize -> {600}, AspectRatio -> 1/GoldenRatio, Frame -> True,
   FrameStyle -> 
    Directive[FontFamily -> "Bookman Old Style", 18, Black], 
   RotateLabel -> False, GridLines -> Automatic};
(*OPCIONS PARTICULARS*)
locals = {FrameLabel -> {{"pH", "" }, {"Vol. valorant [mL]", titol}}, 
   PlotMarkers -> "o"};
ListLinePlot[{data2, Valoracio}, Evaluate[estil], Evaluate[locals]]

Attachments:

exemple.pdf