Off[FindMinimum::"fmmp"];(*???*)
dataPminPmax={{100,625},{80,500},{75,400},{30,300}};(*min. a max. vykony elektraren*)
dataAcka={{500,8,0.004},{400,6.4,0.0048},{600,7.9,0.005},{400,7.5,0.0055}};(*param. nakladovych fci, fix_lin._kvadr.*)
dataCinCout={{3000,1500},{3000,1500},{2000,1000},{1000,800}};(*zapnuti, vypnuti elektraren*)
n=Length[dataAcka];(*pocet elektraren*)
vykony=Table[(Subscript[P, i]),{i,1,n}](*Prirazeni indexu Pi jednotlivym elekrarnam*)
konstantaSkok=(10^6);(*???*)
f[x_]:=(x^2);
pen[x_,a_,b_,kmax_]:=(-kmax((Tanh[1000(x-a)]+1)/2)+kmax((Tanh[1000(x-b)]+1)/2))+kmax(*penalizacni funkce ???*)
ukpen[x_]=pen[x,-20,30,(10^6)];
plPenale=Plot[ukpen[x],{x,-50,50},PlotRange->All,PlotStyle->{Red,Thick},AxesLabel->{"P [MW]","penale"},GridLines->Automatic](*testovaci graf penale*)

penale[i_]:=pen[vykony[[i]],dataPminPmax[[i,1]],dataPminPmax[[i,2]],konstantaSkok];(*ve tvaru napr Pi,intmin,intMax,kmax*)
Plot[penale[1],{(Subscript[P, 1]),0,750},PlotRange->All,PlotStyle->{Red,Thick},AxesLabel->{"P [MW]","penale"},GridLines->Automatic]

cena[i_]:=Module[{Pi},
Pi=vykony[[i]];
Dot[{1,Pi,(Pi^2)},dataAcka[[i]]]];
vahaRovnostDiagram=200;

penaleVykony=Sum[Table[penale[i],{i,1,n,1}][[k]],{k,1,n,1}];(*soucet vsech penale*)
cenaVari=Sum[Table[cena[i],{i,1,n,1}][[k]],{k,1,n,1}];(*soucet vsech elektraren naklady kvadr. fci*)
penaleDodavka[Pdodavka_]:=vahaRovnostDiagram*((((Sum[vykony[[i]],{i,1,n,1}])-Pdodavka))^2)
(*odchylka od denniho diagramu - Secti vsechny vykony, odecti od dodavky, na2 aby to bylo kladne, 200 nejvetsi vaha*)
proFindMinimum[Pdodavka_]:=penaleVykony+penaleDodavka[Pdodavka]+cenaVari

funkcePocPodm[i_]:={vykony[[i]],N[Mean[dataPminPmax[[i]]]]}(*Mean je prumer,N dela desetinne cislo*)
PocPodm=Table[funkcePocPodm[i],{i,1,n,1}]

min[Pdodavka_]:=FindMinimum[proFindMinimum[Pdodavka],PocPodm]
cenaHodina[Pdodavka_]:=cenaVari/.min[Pdodavka][[2]];
datazatizeni={600,800,1100,1400,1000,1200,1000,1400,1600,1600,1400,800,400,200,1400,1200,450,300,400};(*Denni diagram zatizeni po hodinach*)
vysledkyBezKomb=Table[min[datazatizeni[[i]]],{i,1,Length[datazatizeni],1}]//Quiet;

(*optimalizace denniho provozu pres cely den*)

komb={{0,0,0,1},{0,0,1,0},{0,0,1,1},{0,1,0,0},{0,1,0,1},{0,1,1,0},{0,1,1,1},{1,0,0,0},{1,0,0,1},{1,0,1,0},{1,0,1,1},{1,1,0,0},{1,1,0,1},{1,1,1,0},{1,1,1,1}};

nakladyHodina[Pzadany_,kombinace_]:=Module[{penMinMax,penDodavka,pocPodmsKombinacemi,cen,promin,min,prirazeni,cenaCelk},
penMinMax=Dot[Table[penale[i],{i,1,n,1}],kombinace];
penDodavka=(((Dot[vykony,kombinace]-Pzadany))^2);
pocPodmsKombinacemi=Cases[Part[PocPodm,kombinace*Table[i,{i,1,n,1}]],{_,_}];(*vybere dvojice, protoze nulta pozice je list*)
cen=Dot[Table[cena[i],{i,1,n,1}],kombinace];
promin=cen+penDodavka+penMinMax;
min=FindMinimum[promin,pocPodmsKombinacemi]//Quiet;
prirazeni=min[[2]];
cenaCelk=cen/.prirazeni;
{cenaCelk,kombinace,prirazeni,penMinMax+penDodavka/.prirazeni}
];

nakladyHodina[1100,{1,0,1,0}]

vydejHodina[Pzadany_]:=Module[{dats},
dats=Table[nakladyHodina[Pzadany,komb[[i]]],{i,1,Length[komb],1}];
Take[#,3]&/@Sort@Select[dats,#[[4]]<100&](* *)
];

cenaprepnuti[v1_,v2_]:=Sum[If[Sign[(v1-v2)][[i]]<0,dataCinCout[[i,1]],If[Sign[(v1-v2)][[i]]>0,dataCinCout[[i,2]],0]],{i,n}];

start={vydejHodina[datazatizeni[[1]]][[1]],1}

krok[{{c_,komb_,dosP_},n_}]:=Module[{zatizeni,mozne,sprepnutim},
zatizeni=datazatizeni[[n+1]];
mozne=vydejHodina[zatizeni];
sprepnutim=Sort[mozne/.{ci_Real,kombi_List,dosPi_}:>{ci+1*cenaprepnuti[komb,kombi],kombi,dosPi}][[1]];
{sprepnutim,n+1}];

vydejHodina[100];

b+c/.{a->3};

dats=NestList[krok,start,Length[datazatizeni]-1];(*???*)

dataProGraf=Table[vykony/.dats[[i,1,3]],{i,1,Length[datazatizeni]-1}];(*data pro graf*)

BarChart[dataProGraf,ChartLayout->"Stacked",ChartLegends->{"P1","P2","P3","P4"}](*vykresleni grafu*)