Plano de Ensino: Plano de Ensino 2015 Cálculo Numérico
Conteúdo 01: Introdução, erros e Matlab.
Conteúdo sobre EDO´s
//////////////////////////////////////////////////////////////////// /////////Algoritmo para EQUAÇÕES DIFERENCIAIS ORDINÁRIAS //////////////////////////////////////////////////////////////////// %ARQUIVO M chamado 'f.m' function r = f(x,y) %%%Funcoes para Testes dy = 2*x-2; %dy = cos(x); %dy = y/(1+x*x); %%%Exercícios propostos %dy = cos(x); %dy = (x*x)/(1-x); %dy = x*sin(x*x); r = dy;
%ARQUIVO M chamado 'passoRK.m' function r = passoRK(x0,y0,x1) h=x1-x0; k1=h*f(x0,y0); k2=h*f(x0+h/2,y0+k1/2); k3=h*f(x0+h/2,y0+k2/2); k4=h*f(x0+h,y0+k3); y1=y0+(k1/6)+(k2/3)+(k3/3)+(k4/6); r=y1;
%ARQUIVO M chamado 'autoRK.m'
function r = autoRK(x0,y0,x1, precisao)
xm=(x1+x0)/2;
%Tenta passo completo
yt=passoRK(x0,y0,x1);
ym=passoRK(x0,y0,xm);
y1=passoRK(xm,ym,x1);
if(abs(yt-y1)>precisao)
%faltou precisao, dividir em 2 passos e chamar RK recursivamente
ym=autoRK(x0,y0,xm,precisao);
y1=autoRK(xm,ym,x1,precisao);
end
fprintf('RK: x=%.10f y=%.10f\n',x1,y1);
r=y1;
%ARQUIVO M chamado 'nRK.m'
function r = nRK(x0,y0,x1, nPassos)
h=(x1-x0)/nPassos;
x=x0;
y=y0;
sair = false;
i=0;
while not(sair)
xm=x+h;
if xm>x1
xm=x1;
end
ym=passoRK(x,y,xm);
x=xm;
y=ym;
fprintf('RK: x=%.10f y=%.10f\n',x,y);
i=i+1;
if i>=nPassos
sair=true;
end
if x>x1
sair=true;
end
end
r=ym;
%ARQUIVO M chamado 'main.m'
clc;
%%%Funcoes para Testes
x0=0; y0=1; x1=1; % dy = 2*x-2; %Resultado 0.0000
%x0=0; y0=2; x1=pi/2; dy = cos(x); %Resultado 3.0000
%x0=0; y0=1; x1=3.5; % dy = y/(1+(x*x)); %Resultado 2.067 mas algoritmo mostra 3.6418
%%%Exercícios propostos
%x0=0; y0=2; x1=pi/2; % dy = cos(x); %Resultado 3.0000
%x0=2; y0=-1; x1=4; % dy = (x*x)/(1-x); %Resultado -10.0986
%x0=0; y0=0; x1=2*pi; % dy = x*sin(x*x); %Resultado 0.373 mas algoritmo mostra 0.6035
fprintf('Valor inicial:\n');
fprintf('x0 = %.10f\n',x0);
fprintf('y0 = %.10f\n',y0);
fprintf('x1 = %.10f\n',x1);
fprintf('Calculo realizado com autoRK, ou seja, algoritmo aumenta o número de passos automaticamente para chegar na precisão informada:');
precisao=0.00001;
y1 = autoRK(x0,y0,x1,precisao);
fprintf('Resultado = %.10f\n',y1);
numeroDePassos=20;
fprintf('\nCalculo realizado com nRK, ou seja, algoritmo aumenta tem um número de passos definido (atualmente %d passos):',numeroDePassos);
y1 = nRK(x0,y0,x1,numeroDePassos);
fprintf('Resultado = %.10f\n\n',y1);
Conteúdo antigo:
Conteúdo: Método da Dicotomia para encontrar raízes
Baixe a planilha com o método da Dicotomia (Libreoffce) ou Dicotomia (Excel).
Veja também o código em C abaixo:
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
//funcao a ser resolvida
float f(float x)
{
float y;
//y = exp(x)+(x/2.0);
y = pow(x,2.0)-2;
return y;
}
float dicotomia(float xa, float xb, float precisao)
{
float xc, anterior;
anterior = xb;
xc = xa;
while(fabs(xc-anterior)>precisao)
{
anterior = xc;
xc = (xa+xb)/2.0;
if(f(xa)*f(xc)<0.0) xb = xc;
else if(f(xb)*f(xc)<0.0) xa = xc;
}
return xc;
}
int main()
{
float raiz;
raiz = dicotomia(-100, 100, 0.0001);
printf("Raiz = %g\n", raiz);
getchar();//pausa para analisar resultado antes de fechar o programa
return 0;
}
Para finalizar o método da dicotomia, veja o algoritmo abaixo em MATLAB
%ARQUIVO M chamado 'f.m'
function resultado = f(x)
resultado = x^2-2;
%ARQUIVO M chamado 'dicotomia.m'
function r = dicotomia(xa,xb,precisao)
anterior = xb;
xc = xa;
while abs(xc-anterior)>precisao
anterior = xc;
xc = (xa+xb)/2;
if f(xa)*f(xc)<0
xb=xc;
elseif f(xb)*f(xc)<0
xa=xc;
end
end
r = xc;
%ARQUIVO M chamado 'programa.m'
clc
raiz = dicotomia(-100,100,0.0001);
fprintf('Raiz %.4f\n\n',raiz);
%%
%% O resultado com esta função, este intervalo e esta precisão deverá ser:
%% >>Raiz -1.4142
%%
////////////////////////////////////////////////////////////////
Algoritmo de Newton em MATLAB para
encontrar Zeros reais de funções reais
////////////////////////////////////////////////////////////////
clear; clc;
fc = input ('Digite a expressao f(x): '); %Escolha da expressao
epsolon= input ('Digite o valor do criterio de parada (Erro): '); %|xn+1-xn|=<Epsolon
syms x
f = inline(fc);
d=diff(fc, x); %Calcula a derivada
df = inline(d); % cria o objeto inline da derivada
x= input ('Digite a aproximaçao inicial: '); %escolha do valor inicial X0
n=0;
erro = 100;
fprintf('n\t\txi\t\t\tErro\n');
while(erro>=epsolon) %Condiçao para continuar o metodo
fprintf('%d\t\t%.6f\t%.6f\n',n, x, erro)
n=n+1;
x= x - f(x) / df(x);
erro = abs(f(x));
end
fprintf('\n');
%Plotagem do grafico
vi=input('Digite intervalo inicial para o grafico: ');
vf=input('Digite o intervalo final para o grafico: ');
y=[vi : 0.01 : vf];
fplot(f,[vi vf]);
axis([vi vf vi vf]);
grid('on');
zoom(1);
legend (fc);
xlabel ('Eixo X') ;
ylabel ('Eixo Y');
////////////////////////////////////////////////////////////////
Conteúdo sobre Método dos Trapézios
//////////////////////////////////////////////////////////////// Arquivo PDF: Integração Numérica
////////////////////////////////////////////////////////////////
Algoritmo Matlab Método dos Trapézios para Integração Numérica
////////////////////////////////////////////////////////////////
clear; clc;
f = @(x) -x.^2+2*x+2;
syms x
a=0;
b=2;
n=20;
fprintf('Cálculo considerando os seguintes parâmetros: \nIntervalos de a=%.2f até b=%.2f\nNúmero de passos n=%d com a função\n', a, b, n);
f
fprintf('Calculo Analítico');
res = double(int(f(x), 0, 2)) %exibe resposta
fprintf('Calculo pela Regra dos Trapézios Composta Algoritmo 1');
h=(b-a)/n;
res=0;
for k=1:n
x=a+h*k-(h/2);
res=res+f(x)*h;
end
res % Exibe resposta
fprintf('Calculo pela Regra dos Trapézios Composta Algoritmo 2');
h=(b-a)/n;
res=0;
for k=1:(n-1)
x=a+h*k;
res=res+f(x);
end
res=h*((f(a)+f(b))/2)+h*res %exibe resposta
////////////////////////////////////////////////////////////////////////////////////////// Algoritmo Matlab Comparação dos Métodos dos Trapézios e Simpsons para Integração Numérica //////////////////////////////////////////////////////////////////////////////////////////
clear; clc; format long;
f = @(x) -x.^2+2*x+2; syms x; a=0; b=2;
n=10;
fprintf('Cálculo considerando os seguintes parâmetros: \nIntervalos de a=%.2f até b=%.2f\nNúmero de passos n=%d com a função\n', a, b, n);
f
fprintf('Calculo Analítico');
resA = double(int(f(x), a, b)); %exibe resposta
fprintf('\nResultado: %.15f\tErro: %.15f\n\n', resA, abs(resA-resA));
fprintf('Calculo pela Regra dos Trapézios Composta Algoritmo 1');
h=(b-a)/n; res=0;
for k=1:n
x=a+h*k-(h/2);
res=res+f(x)*h;
end
res; % Exibe resposta
fprintf('\nResultado: %.15f\tErro: %.15f\n\n', res, abs(res-resA));
fprintf('Calculo pela Regra dos Trapézios Composta Algoritmo 2');
h=(b-a)/n;
res=0;
for k=1:(n-1)
x=a+h*k;
res=res+f(x);
end
res=h*((f(a)+f(b))/2)+h*res; %exibe resposta
fprintf('\nResultado: %.15f\tErro: %.15f\n\n', res, abs(res-resA));
fprintf('Calculo pela Regra de Simpsons Composta Algoritmo 1');
h=(b-a)/n; res=0;
for k=1:2:(n-1)
x=a+h*k;
res=res+4*f(x);
end
for k=2:2:(n-2)
x=a+h*k;
res=res+2*f(x);
end
res=h*(res+(f(a)+f(b)))/3;
fprintf('\nResultado: %.15f\tErro: %.15f\n\n', res, abs(res-resA));
::::: Método dos Mínimos Quadrados :::::
Segue abaixo o algoritmo em Linguagem C que deve ser convertido para linguagem do MATLAB.
#include
#include
/*Programa REGRESSAO LINEAR
Calcula a reta que melhor se ajusta aos
pontos dados atraves de regressao linear
*/
int N;
float X[30], Y[30];
float A[2][3];
//Ler pontos X e Y do teclado
void LerDados()
{
int i;
printf("Numeros de pontos: n= ");
scanf("%d", &N);
printf("\n\nDigite X e Y separados por espacos.\n\n");
for(i=0; i<N; i++)
{
printf("[%3d] ",i+1);
scanf("%f",&(X[i]));
scanf("%f",&(Y[i]));
}
}
void MontaSistema()
{
int i;
float s;
A[0][1]=0;
A[1][1]=0;
A[0][2]=0;
A[1][2]=0;
for(i=0; i<N; i++)
{
A[0][1] += X[i];
A[1][1] += X[i]*X[i];
A[0][2] += Y[i];
A[1][2] += X[i]*Y[i];
}
A[0][0]=N;
A[1][0]=A[0][1];
}
float det(float a, float b, float c, float d)
{
return (a*d)-(b*c);
}
int main()
{
float a,b;
LerDados();
MontaSistema();
a = det(A[0][2],A[0][1],A[1][2],A[1][1])
/det(A[0][0],A[0][1],A[1][0],A[1][1]);
b = det(A[0][0],A[0][2],A[1][0],A[1][2])
/det(A[0][0],A[0][1],A[1][0],A[1][1]);
printf("g(x)=%10.6f + %10.6f x\n ",a,b);
getch();
return 0;
}
%%%%% SOLUCAO EM MATLAB %%%%%
%%% Arquivo determinante.m
function r = determinante(a,b,c,d)
temp = (a*d)-(b*c);
r = temp;
%%% Arquivo main.m
clc;
%Entrada de dados
%%% Abaixo alimente X e Y e informe em N a quantidade de pontos.
X = [0 1 2 3 4 5 6 7 8 9 ];
Y = [1 3 5 7 9 11 13 15 17 19];
N = 10;
%Monta o sistema
A(1,2)=0;
A(2,2)=0;
A(1,3)=0;
A(2,3)=0;
for i=1:N
A(1,2) = A(1,2)+X(i);
A(2,2) = A(2,2)+X(i)*X(i);
A(1,3) = A(1,3)+Y(i);
A(2,3) = A(2,3)+X(i)*Y(i);
end
A(1,1)=N;
A(2,1)=A(1,2);
%Main
a = determinante(A(1,3),A(1,2),A(2,3),A(2,2))/determinante(A(1,1),A(1,2),A(2,1),A(2,2));
b = determinante(A(1,1),A(1,3),A(2,1),A(2,3))/determinante(A(1,1),A(1,2),A(2,1),A(2,2));
fprintf('g(x)=%10.6f + %10.6f x\n',a,b);
::: Método de GAUSS para resolução de sistemas :::
#include
#include
int N;
float b[30], x[30], a[30][30];
//Algoritmo retirado da obra: Cálculo Numérico, Reinaldo Burian. 2011. Ed. LTC.
//Lê via teclado a ordem do sistema (N) as linhas da matriz de coeficientes e a matriz de resultados
void leMatriz()
{
int i,j;
printf("Ordem do sistema: n= ");
scanf("%d", &N);
printf("\n\nDigite os valores separados por espacos.\n");
for(i=0; i<N; i++)
{
for(j=0; j<N; j++)
scanf("%f", &(a[i][j]));
scanf("%f", &(b[i]));
}
imprimeMatriz();
}
//Método imprimeMatriz() imprime a matriz na tela;
void imprimeMatriz()
{
int i,j;
printf("Ordem do sistema: n= %d\n\n\n", N);
for(i=0; i<N; i++)
{
for(j=0; j<N; j++)
printf(" %10.6f", a[i][j]);
printf(" %10.6f\n", b[i]);
}
printf("\n\n\n ");
getch();
}
//Calcula a matriz das solucoes do sistema ja triangularizado
void retrosubstituicao()
{
int k,j;
float s;
x[N-1]=b[N-1]/a[N-1][N-1];
for(k=N-2; k>=0; k--)
{
s=0;
for(j=k+1; j<N; j++)
s+=a[k][j]*x[j];
x[k]=(b[k]-s)/a[k][k];
}
}
//mostra o resultado da matriz de solucoes
void imprimeSolucao()
{
int i;
printf("\n\nSolucao do sistema linear\n");
for(i=0; i<N; i++)
printf("x(%d)=%10.6g\n",i,x[i]);
printf("\n\n\nProblema encerrado!\n");
}
//Faz a triangularizacao pelo metodo de Gauss
void metodoTriangDeGauss()
{
int i,j,k;
float m;
for(k=0; k<(N-1); k++)
{
for(i=k+1; i<N; i++)
{
m=a[i][k]/a[k][k];//multiplicador
for(j=k+1; j<N; j++)
{
a[i][j] = a[i][j]-m*a[k][j];
imprimeMatriz();
}
b[i]=b[i]-m*b[k];
a[i][k]=0;
imprimeMatriz();
}
}
}
//Abaixo segue o método alternativo para o algoritmo de Gauss
//Calcula o determinante de uma matriz 2x2
float det(float a, float b, float c, float d)
{
return (a*d)-(b*c);
}
//Faz a triangularizacao do sistema pelo metodo de Gauss
void metodoTriangDeGaussAlternativo()
{
int i,j,k;
for(k=0; k<(N-1); k++)
{
for(i=k+1; i<N; i++)
{
for(j=k+1; j<N; j++)
{
a[i][j] = det(a[k][k],a[i][k],a[k][j],a[i][j]);
imprimeMatriz();
}
b[i]=det(a[k][k],a[i][k],b[k],b[i]);
a[i][k]=0;
imprimeMatriz();
}
}
}
int main()
{
// leMatriz();
N=3;
a[0][0]=3; a[0][1]=-2; a[0][2]=5; b[0]=1111;
a[1][0]=6; a[1][1]=-9; a[1][2]=12; b[1]=5111;
a[2][0]=-5; a[2][1]=0; a[2][2]=2; b[2]=1111;
metodoTriangDeGauss();
retrosubstituicao();
imprimeSolucao();
getch();
return 0;
}
//////////////////////////////////
//////////////////////////////////
Correção da última questão da prova:
//////////////////////////////////
//////////////////////////////////
Questão 6 – É possível determinar os juros embutidos em compras a prazo da seguinte forma: verifica-se a cada mês o saldo devedor que no ato da compra do produto equivale ao valor á vista. A cada mês embute-se os juros e subtrai-se a prestação paga. Ao final da operação, passados os meses correspondente ao prazo de pagamento a taxa de juros que faz o montante final ficar nulo é a taxa de juros embutida nas prestações. Sendo pv o valor do produto, pgto o valor da prestação mensal, n o número de prestações e i a taxa de juros, temos a seguinte equação:
Dessa forma, busque a taxa de juros aplicada em uma venda de uma geladeira que á vista custa pv = R$ 2.000,00 que será pago em 24 prestações mensais de pgto = R$ 120,00. A taxa de juros i deve ser buscada com erro menor que 10-6. Utilize um método da dicotomia para resolver este problema e busque justificar os passos (algoritmo) que foi utilizado. Apresente a resolução em planilha e no algoritmo em Matlab.
Resposta: Calculando via Planilha eletronica (Libre Office) temos 0,031748 como resultado correto, se aumentarmos a precisar teremos 0,031491 o que equivale a taxa de juros de 3.15% ao mes.
No Matlab tivemos o resultado Raiz 0.0315 usando a expressão: resultado = (1+x).^24*(2000)-(((1+x).^24-1)/x)*120;
Código em Matlab:
%ARQUIVO M chamado ‘programa.m’
clc
contador = 1;
for contador=1:10
tic % inicia cronometro
raiz = dicotomia(0.01,0.05,0.000001);
fprintf(‘%d: Raiz %.4f’,contador, raiz);
duracao = toc;
fprintf(‘ Duração: %.10f milisegundos\n’, duracao.*1000);
end
%Abaixo segue a impressão do gráfico
x=0.01:0.005:0.04;
%y= (1+x).^84*(25000)-(((1+x).^84-1)/x)*499;
y = (1+x).^24*(2000)-(((1+x).^24-1)/x)*120;
%y=cos(x).*(1/2*exp(x/2));
%cos(x)*1/2*exp(x/2);
%y=x.^3-x.^2+7;
plot(x,y);
%ARQUIVO M chamado ‘f.m’
function resultado = f(x)
%resultado = (1+x).^84*(25000)-(((1+x).^84-1)/x)*499;
resultado = (1+x).^24*(2000)-(((1+x).^24-1)/x)*120;
%resultado = 1+log(x.^3);
%resultado = cos(x)*1/2*exp(x/2);
%resultado = x.^3-x.^2+7;
%ARQUIVO M chamado ‘dicotomia.m’
function r = dicotomia(xa,xb,precisao)
anterior = xb;
xc = xa;
while abs(xc-anterior)>precisao
anterior = xc;
xc = (xa+xb)/2;
if f(xa)*f(xc)<0
xb=xc;
elseif f(xb)*f(xc)<0
xa=xc;
end
end
r = xc;
%%%%%%% CODIGO EM MATLAB PARA O METODO DA TRIANGULARIZACAO DE GAUSS %%%%%%%%%%%%%%%%%%%%%
% Prezados, conforme falamos em sala na aula de 15 de outubro de 2015 segue algoritmo funcional.
% O desafio agora é criar as entradas e saídas de dados personalizadas.
%Método da Triangularização de Gauss
%Atribuição de valores
clear;
clc;
N=3
a = [10 2 1
1 5 1
2 3 10]
b = [7 -8 6]
%Triangularização
for k=1:N-1
for i=k+1:N
m=a(i,k)/a(k,k);%multiplicador
for j=k+1:N
a(i,j)=a(i,j)-m*a(k,j);
a
b
end
b(i)=b(i)-m*b(k);
a(i,k)=0;
a
b
end
end
%Retrosubstituição
k=0; j=0; s=0.0;
x(N)=b(N)/a(N,N);
for k=N-1:-1:1
s=0.0;
for j=k+1:N
s=s+a(k,j)*x(j);
end
x(k) = (b(k)-s)/a(k,k);
end
%Exibição do resultado final
a
b
x
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Exercícios aula dia 29 de outubro:
O algoritmo abaixo em Linguagem C utiliza o método iterativo de Gauss-Seidel para resolução de sistemas
lineares. Converta o Algoritmo para Matlab.
#include
#include
#include
int N;
float b[30], x0[30], x1[30], a[30][30];
float maxDif;
//Algoritmo retirado da obra: Cálculo Numérico, Reinaldo Burian. 2011. Ed. LTC.
//Lê via teclado a ordem do sistema (N) as linhas da matriz de coeficientes e a matriz de resultados
void leMatriz()
{
int i,j;
printf(“Ordem do sistema: n= “);
scanf(“%d”, &N);
printf(“\n\nDigite os valores separados por espacos.\n”);
for(i=0; i<N; i++)
{
for(j=0; j<N; j++)
scanf("%f", &(a[i][j]));
scanf("%f", &(b[i]));
}
}
//Método imprimeMatriz() imprime a matriz na tela;
/*void imprimeMatriz()
{
int i,j;
printf("Ordem do sistema: n= %d\n\n\n", N);
for(i=0; i<N; i++)
{
for(j=0; j<N; j++)
printf(" %10.6f", a[i][j]);
printf(" %10.6f\n", b[i]);
}
printf("\n\n\n “);
getch();
}
*/
/*Calcula a maior diferença entre x1 e x0*/
float maximo()
{
int i;
float d, m=0;
for(i=0; im) m=d;
}
return m;
}
//Mostra resultado da matriz de solucoes
void imprime(float *x, int k, int p)
{
int i;
printf(“\nSolucao do Sistema Linear\n\n”);
for(i=0; i0) printf(“\n\nIteracao %d MaxDif = %g\n\n”,k,maximo());
if(k>0) printf(“\n\nIteracao %d MaxDif = %6.15f\n\n”,k,maximo());
if(p==0) printf(“Problema encerrado.\n”);
printf(“[Enter]”);
getch();
}
/*Faz uma iteração do método*/
void iteracao(int k)
{
int i, j;
for(i=0; i<N; i++) x1[i]=x0[i];
for(i=0; i<N; i++)
{
x1[i] = b[i];
for(j=0; j<N; j++)
if(i!=j)
x1[i] -= a[i][j]*x1[j];
x1[i] /= a[i][i];
}
imprime(x1, k+1, 1);
}
int main()
{
N=3;
a[0][0]=4; a[0][1]=1; a[0][2]=1; b[0]=5;
a[1][0]=3; a[1][1]=1; a[1][2]=6; b[1]=-6.5;
a[2][0]=-2; a[2][1]=5; a[2][2]=1; b[2]=0;
int i, k=0;
//leMatriz();
for(i=0; i<N; i++) x0[i]=0.0;
imprime(x0,k,1);
do
{
iteracao(k);
maxDif = maximo();
for(i=0; i 1.0e-8);
imprime(x1,k,0);
return 0;
}
##################################################################
##################################################################
##################################################################
INTEGRACAO NUMÉRICA
##################################################################
##################################################################
##################################################################
//ALGORTIMO EM C para Trapézios
#include
#include
#include <math.h)
/*FUNCAO A SER INTEGRADA (DEVE SER REESCRITA PARA CADA CASO)*/
float f(float x)
{
float y;
y=1/x;
/*OUTRAS FUNCOES DE EXEMPLO*/
//y=2+2*x-(x*x);
//y=x*sqrt(4-(pow(x,2)));
//y=exp(-pow(x,2)); //exp(x)=e^x (Euler)
return y;
}
/*SOMA DAS AREAS*/
float SomaAreas(float a, float b, float dx)
{
float x, s;
for(x=a+(dx/2), s=0.0; xprecisao);
return s0;
}
int main()
{
float a,b,m;
printf(“\nDigite os limites de integracao:\n\n”);
printf(“a= “);
scanf(“%f”,&a);
printf(“b= “);
scanf(“%f”,&b);
if(a>b)
{
m=a;
a=b;
b=m;
}
printf(“\nResultado: %g\n\n”, Integral(a,b,0.001));
getch();
return 0;
}
##################################################################
##################################################################
##################################################################
Algoritmo em MATLAB Trapézios
%CREDITOS———————————————————————
% AJNeves (DMUA) – Novembro/2005
%—————————————————————————–
%function It=trapezios(f,a,b,N)
%=====================================================
%Integração Numérica
%Regra dos Trapézios Composta
%It=(h/2)[f(a)+2f(x1)+2f(x2)+…+2f(xN-1)+f(b)]
%—————————————————–
%
%Entradas – f função integranda entra como uma
% string ‘f’
% – a é o limite inferior de integração
% – b é o limite superior de integração
% – N é o número de subintervalos
%Saída – It é o valor aproximado do integral
% calculado com a regra dos Trapézios
% composta
%=====================================================
clear; clc;
f = @(x) -x.^2+2*x+2;
syms x
a=0;
b=2;
n=200;
fprintf(‘Cálculo considerando os seguintes parâmetros: \nIntervalos de a=%.2f até b=%.2f\nNúmero de passos n=%d com a função\n’, a, b, n);
f
fprintf(‘Calculo Analítico’);
res = double(int(f(x), 0, 2)) %exibe resposta
fprintf(‘Calculo pela Regra dos Trapézios Composta Algoritmo 1’);
h=(b-a)/n;
res=0;
for k=1:n
x=a+h*k-(h/2);
res=res+f(x)*h;
end
res % Exibe resposta
fprintf(‘Calculo pela Regra dos Trapézios Composta Algoritmo 2’);
h=(b-a)/n;
res=0;
for k=1:(n-1)
x=a+h*k;
res=res+f(x);
end
res=h*((f(a)+f(b))/2)+h*res %exibe resposta
%x=0:0.05:2.5;
%y = @(x) cos(x).*(1/2*exp(x/2))+1;
%area(x,y(x));
%title(‘Gráfico da função f(x)’)
%xlabel(‘x’); ylabel(‘f(x)’); %grid
##################################################################
##################################################################
##################################################################
Algoritmo em MATLAB SIMPSONS
%CREDITOS———————————————————————
% AJNeves (DMUA) – Novembro/2005
%—————————————————————————–
%function Is=simpson(f,a,b,N)
%=================================================
%Integração Numérica
%Regra de Simpson Composta
%Is=(h/3)[f(a)+4f(x1)+2f(x2)+…+4f(xN-1)+f(b)]
%————————————————-
%
%Entradas – f função integranda entra como uma
% string ‘f’
% – a é o limite inferior de integração
% – b é o limite superior de integração
% – N (par) é o número de subintervalos
%Saída – Is é o valor aproximado do integral
% calculado com a regra de Simpson
% composta
%=================================================
clear; clc; format long;
f = @(x) -x.^2+2*x+2; syms x; a=0; b=2;
n=10;
fprintf(‘Cálculo considerando os seguintes parâmetros: \nIntervalos de a=%.2f até b=%.2f\nNúmero de passos n=%d com a função\n’, a, b, n);
f
fprintf(‘Calculo Analítico’);
resA = double(int(f(x), a, b)); %exibe resposta
fprintf(‘\nResultado: %.15f\tErro: %.15f\n\n’, resA, abs(resA-resA));
fprintf(‘Calculo pela Regra dos Trapézios Composta Algoritmo 1’);
h=(b-a)/n; res=0;
for k=1:n
x=a+h*k-(h/2);
res=res+f(x)*h;
end
res; % Exibe resposta
fprintf(‘\nResultado: %.15f\tErro: %.15f\n\n’, res, abs(res-resA));
fprintf(‘Calculo pela Regra dos Trapézios Composta Algoritmo 2’);
h=(b-a)/n;
res=0;
for k=1:(n-1)
x=a+h*k;
res=res+f(x);
end
res=h*((f(a)+f(b))/2)+h*res; %exibe resposta
fprintf(‘\nResultado: %.15f\tErro: %.15f\n\n’, res, abs(res-resA));
fprintf(‘Calculo pela Regra de Simpsons Composta Algoritmo 1’);
h=(b-a)/n; res=0;
for k=1:2:(n-1)
x=a+h*k;
res=res+4*f(x);
end
for k=2:2:(n-2)
x=a+h*k;
res=res+2*f(x);
end
res=h*(res+(f(a)+f(b)))/3;
fprintf(‘\nResultado: %.15f\tErro: %.15f\n\n’, res, abs(res-resA));