File: [local] / OpenXM_contrib / PHC / Ada / Continuation / ts_vlprs.adb (download)
Revision 1.1.1.1 (vendor branch), Sun Oct 29 17:45:22 2000 UTC (23 years, 10 months ago) by maekawa
Branch: PHC, MAIN
CVS Tags: v2, maekawa-ipv6, RELEASE_1_2_3, RELEASE_1_2_2_KNOPPIX_b, RELEASE_1_2_2_KNOPPIX, RELEASE_1_2_2, RELEASE_1_2_1, HEAD Changes since 1.1: +0 -0
lines
Import the second public release of PHCpack.
OKed by Jan Verschelde.
|
with text_io,integer_io; use text_io,integer_io;
with Standard_Floating_Numbers; use Standard_Floating_Numbers;
with Standard_Floating_Numbers_io; use Standard_Floating_Numbers_io;
with Standard_Random_Numbers; use Standard_Random_Numbers;
with Standard_Mathematical_Functions; use Standard_Mathematical_Functions;
with Standard_Integer_Vectors;
with Standard_Integer_Vectors_io; use Standard_Integer_Vectors_io;
with Standard_Floating_Vectors;
with Standard_Floating_Matrices; use Standard_Floating_Matrices;
with vLpRs_Algorithm; use vLpRs_Algorithm;
procedure ts_vlprs is
-- DESCRIPTION :
-- Test on the extrapolation algorithm on power series.
-- GENERATION AND EVALUATION OF A POWER SERIES :
function Power_Series_Exponents
( n : natural ) return Standard_Integer_Vectors.Vector is
-- DESCRIPTION :
-- Asks for a vector of n exponents.
res : Standard_Integer_Vectors.Vector(1..n) := (1..n => 0);
begin
put("Give "); put(n,1); put_line(" integer numbers : ");
for i in res'range loop
get(res(i));
end loop;
return res;
end Power_Series_Exponents;
function Power_Series_Coefficients
( n : natural ) return Standard_Floating_Vectors.Vector is
-- DESCRIPTION :
-- Returns a vector of range 0..n of randomly generated floating-point
-- numbers, which represent the coefficients of a power series.
res : Standard_Floating_Vectors.Vector(1..n);
begin
for i in res'range loop
res(i) := Random;
end loop;
return res;
end Power_Series_Coefficients;
function Eval_Power_Series
( s : double_float; a : Standard_Floating_Vectors.Vector;
m : Standard_Integer_Vectors.Vector ) return double_float is
-- DESCRIPTION :
-- Returns the sum of a_i*s**m(i).
res : double_float := 0.0;
begin
for i in m'range loop
res := res + a(i)*s**m(i);
end loop;
return res;
end Eval_Power_Series;
procedure Eval_Log_Power_Series
( s : in double_float; a : in Standard_Floating_Vectors.Vector;
m : in Standard_Integer_Vectors.Vector;
logs,logx : out double_float ) is
-- DESCRIPTION :
-- On return are the logarithms of s and the result of the power series.
begin
logs := LOG10(s);
logx := LOG10(Eval_Power_Series(s,a,m));
end Eval_Log_Power_Series;
procedure Eval_Log_Power_Series
( s,a : in Standard_Floating_Vectors.Vector;
m : in Standard_Integer_Vectors.Vector;
logs,logx : out Standard_Floating_Vectors.Vector ) is
-- DESCRIPTION :
-- On return are the logarithms of s and the result of the power series.
begin
for i in s'range loop
logs(i) := LOG10(s(i));
logx(i) := LOG10(abs(Eval_Power_Series(s(i),a,m)));
end loop;
end Eval_Log_Power_Series;
-- AUXILIARY :
procedure Write ( l,v : in Standard_Floating_Vectors.Vector;
m : in integer ) is
w,err : double_float;
begin
for i in v'range loop
put(v(i),3,3,3); put(l(i),3,3,3);
w := v(i)/l(i); put(" "); put(w);
err := abs(w - double_float(m));
put(err,3,3,3); new_line;
end loop;
end Write;
-- CALLING the EXTRAPOLATOR :
procedure Extrapolate ( n,m,r : in natural;
a : in Standard_Floating_Vectors.Vector;
exp : in Standard_Integer_Vectors.Vector ) is
s,logs,logx : Standard_Floating_Vectors.Vector(0..m-1);
srp,dsp : Standard_Floating_Vectors.Vector(1..r-1) := (1..r-1 => 0.0);
p : Standard_Floating_Vectors.Vector(0..r-1) := (0..r-1 => 0.0);
l,v : Standard_Floating_Vectors.Vector(0..r) := (0..r => 0.0);
rt1,rt2 : Matrix(1..r-1,1..r-1);
begin
put("Give an decreasing sequence of "); put(m,1);
put_line(" floating-point numbers :");
for i in s'range loop
get(s(i));
end loop;
Eval_Log_Power_Series(s,a,exp,logs,logx);
put_line("The logs of the s-values : ");
for i in logs'range loop
put(logs(i)); new_line;
end loop;
put_line("The logs of the x-values : ");
for i in logx'range loop
put(logx(i)); new_line;
end loop;
rt1(1,1) := 0.0; rt2(1,1) := 0.0;
vlprs_full(r,s,logs,logx,srp,dsp,p,l,v,rt1,rt2);
put_line("The extrapolated values : ");
Write(l,v,exp(exp'first));
put_line("The errors of the extrapolated values : ");
vlprs_pipe(Standard_Output,r,s,logs,logx,srp,dsp,p,l,v,rt1,rt2);
put_line("The extrapolated values : ");
Write(l,v,exp(exp'first));
end Extrapolate;
procedure Test_vlprs is
n,m,r : natural;
begin
put("Give the order of the extrapolator : "); get(r);
loop
put("Give the number of points : "); get(m);
exit when m > r;
put("should be larger than "); put(r,1); put_line(". Please retry.");
end loop;
put("Give the length of the power series : "); get(n);
declare
exp : Standard_Integer_Vectors.Vector(1..n) := Power_Series_Exponents(n);
pcf : Standard_Floating_Vectors.Vector(1..n)
:= Power_Series_Coefficients(n);
begin
put("the exponents : "); put(exp); new_line;
put_line("the coefficients : ");
for i in pcf'range loop
put(pcf(i)); new_line;
end loop;
Extrapolate(n,m,r,pcf,exp);
end;
end Test_vlprs;
begin
Test_vlprs;
end ts_vlprs;