Integração Numérica - 1ª regra de Simpson

O programa calcula a integral definida (de a até b) de uma função contínua neste intervalo (Integral[a,b]f(x)dx). O usuário deve digitar no código fonte a equação y = f(x) que deseja integrar, além dos limites da integração em x, e depois compilar em qualquer compilador pascal para ver o resultado.

Exemplo: para integrar a função y = x^3 - e^x, no intervalo de x = -3 até 5 o usuário irá digitar antes de compilar:

...

y := x*x*x - exp(x);

...

lim_inferior := -3;

lim_superior := 5;

...

Observação: Quanto maior o valor de n (número de divisões do intervalo de integração) maior será a exatidão do resultado da integração numérica, mas, em contrapartida, mais cálculos serão feitos, tornando o processo mais lento.

O local onde o usúario deve entrar com os dados esta destacado em AZUL:


program integracao_1_Simpson;
uses crt;

function y (x:real):real;
begin
y := x*x*x - exp(x); //digitar a funcao aqui//
end;

var
lim_inferior, lim_superior, h, n, dx, xj, xjm, cont, I : real;

begin

lim_inferior := -3; //digitar o limite inferior da integracao//
lim_superior := 5; //digitar o limite superior da integracao//

//ajusta algumas variaveis//
n := 1000; //numero de divisoes do intervalo de integracao - n deve ser um inteiro positivo//
dx := (lim_superior - lim_inferior)/n;
h := dx/2;
I := y(lim_inferior) + 4*y(lim_superior-h) + y(lim_superior);
cont := 1;

if n > 1 then
begin

xj := lim_inferior + h;
xjm := lim_inferior + dx;

//integracao numerica//

repeat

I := I + 4*y(xj) + 2*y(xjm);
cont := cont + 1;

xj := xj + dx;
xjm := xjm + dx;

until cont = n;

end;

I := (h/3)*I;

//resultados//
writeln ('resultado da integracao = ', I:0:10);
writeln;
writeln ('Digite uma tecla para sair');
readkey;
end.

Integração Numérica - 1ª regra de Simpson - Cálculo da Área, Centróide e Volume de Revolução de uma área gerada a partir de uma equação.

O programa abaixo calcula a área abaixo da curva de uma equação qualquer, além do centróide dessa área e o volume de revolução que essa área gera a partir de uma revolução completa em torno do eixo x e do eixo y.

O usuário deve digitar no código fonte a equação em função de x e em função de y, além dos limites da integração em x e y (valores máximos e mínimos de x e y do intervalo escolhido) e depois compilar em qualquer compilador pascal para ver os resultados.

Exemplo: a função y = 16 - 4x, no intervalo de x = 0 até 4. O usuário irá digitar antes de compilar:

...

funcao_x := 16 - 4*x; //y = 16 - 4x//

funcao_y := (16 - y)/4; //x = (16 - y)/4//

...

x_min := 0; 
x_max := 4; 
...

Observação: Quanto maior o valor de n (número de subdivisões do intervalo de integração) maior será a exatidão do resultado da integração numérica, mas, em contrapartida, mais cálculos serão feitos, tornando o processo mais lento.

O Local onde o usúario deve entrar com os dados esta destacado em AZUL:


program centroide_1_Simpson;
uses crt;

function funcao_x (x:real):real;
begin
funcao_x := 16 - 4*x; //y = f(x). digitar aqui a equacao em funcao de x//
end;

function funcao_y (y:real):real;
begin
funcao_y := (16 - y)/4 ; //x = f(y). digitar aqui a equacao em funcao de y//
end;

var
aux, x_min, x_max, y_min, y_max, dx, dy, xj, xjm, yj, yjm, S, ds, area, cont,
momento_x, momento_y, c_x, c_y, vol_x, vol_y, area_x, area_y, n, h, z, sinal, divisoes : real;
begin

textbackground (black);
textcolor (white);
clrscr;

//limites da figura (limites da integracao) no eixo x - minimo e maximo//

//entrar com esses valores//
x_min := 0;  //x minimo//
x_max := 4; //x maximo//


y_min := funcao_x(x_min);
y_max := funcao_x(x_max);
if Abs(y_min) >= Abs(y_max) then
begin
aux := y_max;
y_max := y_min;
y_min := aux;
z := 0;
sinal := -1;
end
else
begin
z := x_max;
sinal := 1;
end;

//calculo da area da figura//

//ajusta algumas variaveis//
n := 10000; //n tem que ser par//
dx := (x_max - x_min)/n;
h := dx/2;
xj := x_min;
xjm := x_min + dx;
yj := (xj + h);
area := funcao_x(xj) + 4*funcao_x(yj) + funcao_x(xjm);
cont := 1;

//integracao numerica - 1 regra de Simpson//

repeat

xj := xjm;
yj := xjm + h;
xjm := xjm + dx;

area := area + Abs(funcao_x(xj) + 4*funcao_x(yj) + funcao_x(xjm));
cont := cont + 1;

until cont = n;

area := (h/3)*area;


//momento estatico (primeira ordem) em relacao a y//

//ajusta algumas variaveis//
n := 10000; //n tem que ser par//
dx := (x_max - x_min)/n;
h := dx/2;
xj := x_min;
xjm := x_min + dx;
yj := (xj + h);
S := funcao_x(xj)*xj + 4*funcao_x(yj)*yj + funcao_x(xjm)*xjm;
cont := 1;

//integracao numerica - 1 regra de Simpson//

repeat

xj := xjm;
yj := xjm + h;
xjm := xjm + dx;

S := S + funcao_x(xj)*xj + 4*funcao_x(yj)*yj + funcao_x(xjm)*xjm;
cont := cont + 1;

until cont = n;

momento_y := (h/3)*S;


//momento estatico (primeira ordem) em relacao a x//

//ajusta algumas variaveis//
n := 10000; //n tem que ser par//
dy := (y_max - y_min)/n;
h := dy/2;
yj := y_min;
yjm := y_min + dy;
xj := (yj + h);
S := ((z - funcao_y(yj))*yj + 4*(z - funcao_y(xj))*xj + (z - funcao_y(yjm))*yjm)*sinal;
cont := 1;

//integracao numerica - 1 regra de Simpson//

repeat

yj := yjm;
xj := yjm + h;
yjm := yjm + dy;

S := S + ((z - funcao_y(yj))*yj + 4*(z - funcao_y(xj))*xj + (z - funcao_y(yjm))*yjm)*sinal;
cont := cont + 1;

until cont = n;

momento_x := (h/3)*S;


//coordenadas do centroide//
c_x := momento_y / area;
c_y := momento_x / area;

//volume de revolucao//
vol_x := Abs(2*pi*c_y*area);
vol_y := Abs(2*pi*c_x*area);

//area de revolucao//

//ajusta algumas variaveis//
divisoes := 1000000;
dx := (x_max - x_min)/divisoes;
xjm := x_min;
S := 0;
cont := 0; //contador//

//comprimento da curva//
repeat

xj := xjm;
xjm := xjm + dx;

yj := funcao_x(xj);
yjm := funcao_x(xjm);

ds := sqrt(sqr(xjm - xj) + sqr(yjm - yj)); //trecho da curva//

S := S + ds; //somatorio dos trechos da curva//

cont := cont + 1;

until cont = divisoes;

area_x := Abs(2*pi*c_y*S);
area_y := Abs(2*pi*c_x*S);


//resultados//

writeln ('area = ', area:0:4);
writeln;
writeln ('momento estatico em x = ', momento_x:0:4);
writeln ('momento estatico em y = ', momento_y:0:4);
writeln;
writeln ('coordenada x do centroide = ', c_x:0:4);
writeln ('coordenada y do centroide = ', c_y:0:4);
writeln;
writeln ('area da superficie de revolucao em torno do eixo x = ', area_x:0:4);
writeln ('volume do corpo de revolucao em torno do eixo x = ', vol_x:0:4);
writeln ('(para um revolucao completa (de 360 graus = 2pi) em torno do eixo x)');
writeln;
writeln ('area da superficie de revolucao em torno do eixo y = ', area_y:0:4);
writeln ('volume do corpo de revolucao em torno do eixo y = ', vol_y:0:4);
writeln ('(para um revolucao completa (de 360 graus = 2pi) em torno do eixo y)');
writeln;
writeln ('equacoes base para o centroide:');
writeln;
writeln ('coordenada x do centroide = ', (c_x/x_max):0:4, ' * x maximo');
writeln ('coordenada y do centroide = ', (c_y/y_max):0:4, ' * y maximo');
writeln;
writeln ('Digite uma tecla para sair');
readkey;
end.