I'm trying to implement the rectangle method in Pascal. The point is, I'm getting wrong results. I'm getting 0, while the same code in C gives me good results. Why is that? Thanks.
Pascal:
Program HelloWorld(output);
function degrees2radians(x: real) : real;
var
result: real;
begin
result := x * 3.14159 / 180.0;
end;
function fun1(x: real) : real;
var
result: real;
begin
x := degrees2radians(x);
result := Sin(x);
end;
function rect_method(a: real; b: real; n: integer) : real;
var
result: real;
h: real;
integral: real;
i : integer;
begin
integral := 0;
i := 0;
h := (a - b) / n;
for i := 1 to n do
begin
integral := integral fun1(b i*h)*h;
end;
result := integral;
end;
var
result: real;
begin
result := rect_method(1.0, -2.0, 3);
writeln('result = ', result);
end.
And C (which works: https://onlinegdb.com/ubuNQInB2):
#include <iostream>
#include <cstdlib>
#include <cmath>
double degrees2radians(double x)
{
return x * 3.14159 / 180;
}
double f1(double x)
{
x = degrees2radians(x);
return sin(x);
}
double rectangle_method(double a, double b, int n)
{
float h, integral = 0;
h = (a - b) / (double) n;
for (int i=1; i<=n; i )
integral = f1(b i*h)*h;
return integral;
}
int main()
{
std::cout << rectangle_method(1, -2, 3) << "\n";
return 0;
}
CodePudding user response:
In your Pascal code result
is a local variable, which has nothing to do with the special identifier called result
that you want to use:
function degrees2radians(x: real) : real;
var
result: real;
begin
result := x * 3.14159 / 180.0;
end;
You should remove result: real;
in all your functions.
However, "fixed" code produces zero anyway: https://onlinegdb.com/_M2pGk134, which is actually correct. That is, this code should return zero.
I rewrote this code in Julia (the language doesn't matter, the point is that it supports high-precision floating-point types out of the box; BigFloat
is such a high-precision float), and the result is still zero:
degrees2radians(x::Real)::Real = x * 3.14159 / 180
f1(x::Real) = sin(degrees2radians(x))
function rectangle_method(a::Real, b::Real, n::Integer)
integral = 0
h = (a - b) / n;
for i in 1:n
@show i f1(b i*h)*h
integral = f1(b i*h)*h;
end
return integral;
end
@show rectangle_method(BigFloat("1"), -BigFloat("2"), 3)
The output looks like this:
i = 1
f1(b i * h) * h = -0.01745239169736329549313317881927049082689975241899824937795751704833664190229729
i = 2
f1(b i * h) * h = 0.0
i = 3
f1(b i * h) * h = 0.01745239169736329549313317881927049082689975241899824937795751704833664190229729
rectangle_method(BigFloat("1"), -(BigFloat("2")), 3) = 0.0
So, f1(b 2 * h) * h
is zero, and f1(b 1 * h) * h
is exactly (look at the amount of digits after the decimal point!) the negative of f1(b 3 * h) * h
, so they cancel out in the integral
sum, resulting in zero.