f[x_]:=x^2+2*x*Sin[x]-Cos[x]-1;
n=5;
a=0;
b=π;
ya=0;
yb=0;
a0=-1;
p[x_]:=x;
q[x_]:=x;
presneriesenie[x_]:=Sin[x]+x;

Okrsamoadjnp[n_,xleft_,xright_,Leftbv_,Rightbv_,p_,q_,f_,pres_]:=
(
h=N[(xright-xleft)/(n+1)];
ps=Join[{f[xleft+h]*h^2+p[xleft+h/2]*Leftbv},Table[h^2*f[xleft+i*h],{i,2,n}],{f[xleft+(n+1)*h]*h^2+p[xright-h/2]+Rightbv*2*h}];
diag=Table[p[xleft+i*h-h/2]+p[xleft+i*h+h/2]+h^2*q[xleft+i*h],{i,1,n+1}];
hdiag=Table[-p[xleft+i*h+h/2],{i,1,n}];
ddiag=Join[Table[-p[xleft+i*h+h/2],{i,1,n-1}],{-(p[xleft+(n+1)*h-h/2]+p[xleft+(n+1)*h+h/2])}];
A=SparseArray[{Band[{1,2}]->hdiag,Band[{1,1}]->diag,Band[{2,1}]->hdiag},{n,n}];
riesenie=LinearSolve[A,N[ps]];
xs=Table[xleft+(i-1)*h,{i,1,n+2}];
ys=Join[{Leftbv},Table[riesenie[[i]],ys[[i]]],{i,1,n+2}];
graf=Table[{xs[[i]],ys[[i]]},{i,1,n+2}];vysledky=Table[{xs[[i]],ys[[i]],pres[xs[[i]]],Abs[(ys[[i]]-pres[xs[[i]]])]},{i,1,n+2}];
Print[MatrixForm[A]];
riesenie=LinearSolve[A,N[ps]];
)

Okrsamoadjnp[n,a,b,ya,yb,p,q,f,presneriesenie]