I : integer;

begin

x2:=x*x;

sum:=x;

term:=x;

i:=0;

repeat

i:=i+1;

sum1:=sum;

term:=2.0*term*x2/(1.0+2.0*i);

sum:=term+sum1

until term<tol*sum;

erf:=2.0*sum*exp(-x2)/sqrtpi

end; { erf }

function erfc(x: real): real;

{ complement of error function }

const sqrtpi = 1.7724538;

terms = 12;

var x2,u,v,sum : real;

I : integer;

begin

x2:=x*x;

v:=1.0/(2.0*x2);

u:=1.0+v*(terms+1.0);

for i:=terms downto 1 do

begin

sum:=1.0+i*v/u;

u:=sum

end;

erfc:=exp(-x2)/(x*sum*sqrtpi)

end; { ercf }

begin { main }

ClrScr;

done:=false;

writeln;

repeat

write('Arg? ');

readln(x);

if x<0.0 then done:=true

else

begin

if x=0.0 then

begin

er:=0.0;

ec:=1.0

end

else

begin

if x<1.5 then

begin

er:=erf(x);

ec:=1.0-er

end

else

begin

ec:=erfc(x);

er:=1.0-ec

end { if }

end;

writeln('X= ',x:8:4,' Erf= ',er:12:8,', Erfc= ',ec:12)

end { if }

until done

end.

Программа 25. Решение системы уравнений методом Гаусса-Зейделя.

program gausid;

{ pascal program to perform simultaneous solution }

{ by Gauss-Seidel }

{ procedure SEID is included }

const maxr = 8;

maxc = 8;

type ary = array[1..maxr] of real;

arys = array[1..maxc] of real;

ary2s = array[1..maxr,1..maxc] of real;

var y : ary;

coef : arys;

a : ary2s;

n,m : integer;

first,

error : boolean;

procedure get_data

(var a : ary2s;

var y : ary;

var n,m: integer);

{ get values for n and arrays a,y }

var i,j : integer;

begin

writeln;

repeat

write('How many equations? ');

readln(n);

if first then first:=false else ClrScr

until n<maxr;

m:=n;

if n>1 then

begin

for i:=1 to n do

begin

writeln('Equation',i:3);

for j:=1 to n do

begin

write(j:3,':');

read(a[i,j])

end;

write(' C:');

read(y[i]);

readln { clear the line }

end;

writeln;

for i:=1 to n do

begin

for j:=1 to m do

write(a[i,j]:7:4,' ');

writeln(':',y[i]:7:4)

end;

writeln

end { if n>1 }

else if n<0 then n:=-n;

m:=n

end; { procedure get_data }

procedure write_data;

{ print out the answers }

Var I : integer;

begin

for i:=1 to m do

write(coef[i]:9:5);

writeln

end; { write_data }

procedure seid

(a : ary2s;

y : ary;

var coef: arys;

ncol : integer;

var error: boolean);

{ matrix solution by Gauss Seidel }

const tol = 1.0E-4;

max = 100;

var done : boolean;

i,j,k,l,n: integer;

nextc,hold,

sum,lambda,

ab,big : real;

begin

repeat

write('Relaxation factor? ');

readln(lambda)

until (lambda<2) and (lambda>0.0);

error:=false;

n:=ncol;

for i:=1 to n-1 do

begin

big:=abs(a[i,i]);

l:=i;

for j:=i+1 to n do

begin

{ search for largest element }

ab:=abs(a[j,i]);

if ab>big then

begin

big:=ab;

l:=j

end

end; { j-loop }

if big=0.0 then error:=true

else

begin

if l<>i then

begin

{ interchange rows to put }

{ largest element on diagonal }

for j:=1 to n do

begin

hold:=a[l,j];

a[l,j]:=a[i,j];

a[i,j]:=hold

end;

hold:=y[l];

y[l]:=y[i];

y[i]:=hold

end { if l<>i }

end { if big }

end; { i-loop }

if a[n,n]=0.0 then error:=true

else

begin

for i:=1 to n do

coef[i]:=0.0; { initial guess }

i:=0;

repeat

i:=i+1;

done:=true;

for j:=1 to n do

begin

sum:=y[j];

for k:=1 to n do

if j<>k then

sum:=sum-a[j,k]*coef[k];

nextc:=sum/a[j,j];

if abs(nextc-coef[j])>tol then

begin

done:=false;

if nextc*coef[j]<0.0 then

nextc:=(coef[j]+nextc)*0.5

end;

coef[j]:=lambda*nextc+(1.0-lambda)*coef[j];

writeln(i:4,',coef(',j,')=',coef[j])

end { j-loop }

until done or (i>max)

end; { if a[n,n]=0 }

if i>max then error:=true;

if error then writeln('ERROR: Matrix is singular')

end; { SEID }

begin { MAIN program }

first:=true;

ClrScr;

writeln;

writeln('Simultaneous solution by Gauss-Seidel');

repeat

get_data(a,y,n,m);

if n>1 then

begin

seid(a,y,coef,n,error);

if not error then write_data

end

until n<2

end.

Программа 26. Параболическая интерполяция с помощью МНК.

program least1;

{ Pascal Program to perform a liner least-squares fit using a parabolic }

{ curve. Separate procedure PLOT needed }

const maxr = 20;

maxc = 3;

type ary = array[1..maxr] of real;

arys = array[1..maxc] of real;

ary2s = array[1..maxc,1..maxc] of real;

var x,y,y_calc : ary;

coef : arys;

nrow,ncol : integer;

correl_coef : real;

procedure get_data(var x,y: ary;

var nrow: integer);

{ get values for n and arrays x,y }