Interb-04-03-21, R

Interb-04-03-21

Language: Layout:

f<-function(s){log(s^2*sin(s) + cos(s))}
curve(f,0,1)

x=seq(0,1,by=0.1);x
y=f(x);y

plot(x,y)  # Supomos y=f(x) e que f(x_i)=p(x_i), em que p(s)=a_0+a_1s+a_2s^2+a_3s^3+a_4s^4 é um polinômio.

# Podemos resolver esse problema resolvendo o sistema linear 
# Au=y
# em que a matriz A tem coluna i dada pelo vetor (1,x_i,x_i^2,...,x_i^(n-1)).

n=length(x)
A=matrix(0,n,n)

for ( i in 1:n){
for (j in 1:n){
A[i,j]=x[i]^(j-1)
}
}

gradconj<-function(A,b,m){      # Início do método
B=t(A)%*%A; bb=t(A)%*%b           # Adaptando o sistema para A*Au=A*v, ou Bu=b, em que B=A*A é positiva definida e b=A*v.
u=0*b
w=B%*%u-bb                        # w0
v=w                              # d1
alpha= t(w)%*%w /(t(B%*%w) %*%w); alpha=alpha[1]     # Escolha de h0
u=u-alpha*w                       # u1

for ( j in 1:m){
w=(B%*%u-bb)
beta=t(B%*%w) %*%v/(t(B%*%v) %*%v); beta=beta[1]
v=w-beta*v                   # Método de Euler
alpha= t(v) %*%w/(t(B%*%v) %*%v);alpha=alpha[1]
u=u-alpha*v                      # Método quase Euler (integral)
}
return(u)
}

print(" Aproximação para u");u=gradconj(A,y,50);u

print("Erro em Au-y"); A%*%u-y
#---------------------------------------------------------------------------------------

poli<-function(s){
q=u[1]
for (i in 2:n){q=q+u[i]*s^(i-1)}
q
}

curve(poli,min(x),max(x))                           # Teste de ajuste
points(x,y,col="red")

Erro<-function(s){f(s)-poli(s)}
curve(Erro,0,1)

[ + ] Show input
Absolute running time: 0.94 sec, cpu time: 0.95 sec, memory peak: 51 Mb, absolute service time: 1,14 sec

[1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 [1] 0.000000000 -0.004005512 -0.012059068 -0.018231888 -0.016771937 [6] -0.002564339 0.028205368 0.077432091 0.144805383 0.228015476 [11] 0.323367668 [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [1,] 1 0.0 0.00 0.000 0.0000 0.00000 0.000000 0.0000000 0.00000000 [2,] 1 0.1 0.01 0.001 0.0001 0.00001 0.000001 0.0000001 0.00000001 [3,] 1 0.2 0.04 0.008 0.0016 0.00032 0.000064 0.0000128 0.00000256 [4,] 1 0.3 0.09 0.027 0.0081 0.00243 0.000729 0.0002187 0.00006561 [5,] 1 0.4 0.16 0.064 0.0256 0.01024 0.004096 0.0016384 0.00065536 [6,] 1 0.5 0.25 0.125 0.0625 0.03125 0.015625 0.0078125 0.00390625 [7,] 1 0.6 0.36 0.216 0.1296 0.07776 0.046656 0.0279936 0.01679616 [8,] 1 0.7 0.49 0.343 0.2401 0.16807 0.117649 0.0823543 0.05764801 [9,] 1 0.8 0.64 0.512 0.4096 0.32768 0.262144 0.2097152 0.16777216 [10,] 1 0.9 0.81 0.729 0.6561 0.59049 0.531441 0.4782969 0.43046721 [11,] 1 1.0 1.00 1.000 1.0000 1.00000 1.000000 1.0000000 1.00000000 [,10] [,11] [1,] 0.000000000 0.0000000000 [2,] 0.000000001 0.0000000001 [3,] 0.000000512 0.0000001024 [4,] 0.000019683 0.0000059049 [5,] 0.000262144 0.0001048576 [6,] 0.001953125 0.0009765625 [7,] 0.010077696 0.0060466176 [8,] 0.040353607 0.0282475249 [9,] 0.134217728 0.1073741824 [10,] 0.387420489 0.3486784401 [11,] 1.000000000 1.0000000000 [1] " Aproximação para u" [,1] [1,] 1.295061e-08 [2,] -6.884064e-05 [3,] -4.984836e-01 [4,] 9.871671e-01 [5,] -2.456162e-02 [6,] 1.590358e-01 [7,] -1.578435e-01 [8,] -3.694218e-01 [9,] -4.243405e-03 [10,] 4.128681e-01 [11,] -1.810805e-01 [1] "Erro em Au-y" [,1] [1,] 1.295061e-08 [2,] -8.869453e-08 [3,] 2.485736e-07 [4,] -3.407009e-07 [5,] 1.507262e-07 [6,] 2.365295e-07 [7,] -4.657695e-07 [8,] 3.833775e-07 [9,] -1.768744e-07 [10,] 4.475192e-08 [11,] -4.869432e-09