kert:=proc(xv)
	local  n,j,b,ma ;
	n:=Dimension(xv);
	b:=Vector(n);
	for j to n do 
		b[j] := 1/j; 
	end do;
	ma := Transpose(VandermondeMatrix(xv));
	LinearSolve(ma, b);
end proc;


inteapp:=proc(fu,ke,x)
	local n,va,j ;
	n:=Dimension(xv);
	va:=0;
	for j to n do 
		va := va+ke[j]*fu(x[j]); 
	end do;
	va;
end proc;

inteapp1:=proc(fu,x)
	local n,va,j,ke ;
	n:=Dimension(xv);
	ke:=kert(x);
	va:=0;
	for j to n do 
		va := va+ke[j]*fu(x[j]); 
	end do;
	va;
end proc;


kertnc:=proc(n)
	local  j,b,x,ma ;
	b:=Vector(n+1);
	x:=Vector(n+1);
	for j to n+1 do 
		b[j] := 1/j; 
		x[j]:=evalf((j-1)/n);
	end do;
	ma := Transpose(VandermondeMatrix(x));
	LinearSolve(ma, b);
end proc;

kertpistga:=proc(n)
	local po,pist  ;
	po:=subs([s=2*x-1],P(n,s));
	pist:=fsolve(po);
	pist:=convert([pist], Vector);
	[kert(pist),pist];
end proc;


intbe:=proc(fu,n)
	local va,j;
	va:=0;
	for j to n+1 do 
		va:=va+fu((j-1)*h);
	end do;
	va:=evalf(va/(n+1));
end proc;

intnc:=proc(fu,n)
	local ke,va,h,j;
	ke:=kertnc(n);
	va:=0;
	h:=evalf(1/n);
	for j to n+1 do 
		va:=va+ke[j]*fu((j-1)*h);
	end do;
	va;
end proc;
	

intga:=proc(fu,n)
	local kepi, ke, pis ,va,h,j;
	kepi:=kertpistga(n);
	ke:=kepi[1];
	pis:=kepi[2];
	va:=0;
	for j to n do 
		va:=va+ke[j]*fu(pis[j]);
	end do;
	va;
end proc;
	
	
inttrap:=proc(fu,n)
	local ke,va,h,j;
	va:=evalf(fu(0)+fu(1));
	h:=evalf(1/n);
	for j from 2 to n do 
		va:=va+2*evalf(fu((j-1)*h));
	end do;
	va:=h*va/2;
end proc;


intsimp:=proc(fu,n)
	local ke,va,h,j,m;
	va:=evalf(fu(0)+fu(1));
	h:=evalf(1/n);
	m:=n/2;
	for j from 1 to m-1 do 
		va:=va+4*fu((2*j-1)*h);
		va:=va+2*fu((2*j)*h);
	end do;
	va:=va+4*fu((2*m-1)*h);
	va:=h*va/3;
end proc;


inttrap2d:=proc(fu,n)
	local i,ke,va,h,j,m,mc;
	h:=evalf(1/n);
	ke := Vector(n+1); 
	ke[1] := 1; 
	ke[n+1] := 1; 
	for j from 2 to n do 
		ke[j] := 2 
	end do; 
	mc := OuterProductMatrix(ke, ke);
	
	va := 0; 
	for i to n+1 do 
		for j to n+1 do 
			va := va+evalf(mc[i, j]*fu((i-1)*h, (j-1)*h)) 
		end do;
	end do; 
	va := h^2*va/4;
end proc;


intsimp2d:=proc(fu,n)
	local i,ke,va,h,j,m,mc;
	m:=n/2;
	h:=evalf(1/n);
	ke := Vector(n+1); 
	ke[1] := 1; 
	ke[n+1] := 1; 
	for j to m do 
		ke[2*j] := 4 
	end do; 
	for j to m-1 do 
		ke[2*j+1] := 2 
	end do; 
	mc := OuterProductMatrix(ke, ke);
	
	va := 0; 
	for i to n+1 do 
		for j to n+1 do 
			va := va+evalf(mc[i, j]*fu((i-1)*h, (j-1)*h)) 
		end do;
	end do; 
	va := h^2*va/9;
end proc;

	
virhetr:=proc(fu,m,ta)
	local  k,vi,ap,n;
	vi:=Matrix(m+1,2);
	for k from 0 to m do  
		n:= 4*2^k;
		ap:=inttrap(fu,n);
		vi[k+1,1]:=evalf(1/n);
		vi[k+1,2]:=abs((ta-ap)/ta);
	end do;
	vi;
end proc;

virhesi:=proc(fu,m,ta)
	local  k,vi,ap,n;
	vi:=Matrix(m+1,2);
	for k from 0 to m do  
		n:= 4*2^k;
		ap:=intsimp(fu,n);
		vi[k+1,1]:=evalf(1/n);
		vi[k+1,2]:=abs((ta-ap)/ta);
	end do;
	vi;
end proc;

virhega:=proc(fu,m,ta)
	local  k,vi,ap,n;
	vi:=Matrix(m+1,2);
	for k from 0 to m do  
		n:= 4*2^k;
		ap:=intga(fu,n);
		vi[k+1,1]:=evalf(1/(n+1));
		vi[k+1,2]:=abs((ta-ap)/ta);
	end do;
	vi;
end proc;