F<-function(s,u){
    1+u^2
}

# Condição inicial y(0)=0, y'=1+y^2=F(t,y). 
n=10
a=0; b=1
h=(b-a)/n 

t=seq(a,b,by=h)
w=0*t 

w[1]=0 # condição incial y(0)=0

for (i in 1:n){ # Runge-Kutta4
k1=F(t[i],w[i])
k2=F(t[i]+h/3,w[i]+h*k1/3)
k3=F(t[i]+2*h/3,w[i]-h*k1/3+h*k2)
k4=F(t[i]+h,w[i]+h*(k1-k2+k3))
 w[i+1]=w[i]+h*(k1+3*(k2+k3)+k4)/8
}

# Sabemos que y(t)=tan(t)
f<-function(s){tan(s)}
curve(f,a,b,col="blue",ylab="tg(x)")

points(t,w,col="red")

#------------ Interpolação 


n=length(w)

L<-function(i,s){
p=1
for (j in 1:n){
if ( i!=j){p=p*(s-t[j])/(t[i]-t[j])}
}
p
}

poliLag<-function(s){
p=0
for ( i in 1:n){p=p+w[i]*L(i,s)}
p
}

curve(poliLag,a,b)
points(t,w,col="red")

Erro<-function(s){f(s)-poliLag(s)}
curve(Erro,a,b,ylab="tg(x)-poliLag(x)")
points(t,f(t)-w,col="red")