Технический форум - Нужна помощь с программой расчета колебаний струны

Доброго времени суток! Мне срочно нужна помощь. Я пишу программу. расчитывающую колебания струны. на выходе файл с таблицей. ось узлов и ось времени. Я мат. расчеты, в их точности я уверен. Ниже код, в нем не очень трудно разобраться. 2 массива - для таблицы по времени и узлам. формула расчета верна *(((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + ( (f *(Ht * Ht)) + (2 * Arr[i,j] - Arr[i-1,j]))* могу и фотографии расчетов отправить. я постарался расписать действия в программе. но значения не адекватные выдаёт. в общем Очень прошу помощи. завтра буквально нужно её представить.

program struna;
type Mass = Array[0..1000, 0..1000] Of Real;
label M1;
var
S : Text;
f: Real;
Q : real;
Fi : real;
d : real;
Ht: Real;
Hx: Real;
ft: Real;
a: Real;
t: integer;
i, j : integer;
n, m: integer;
Arr: Mass;
fCSV: TEXT;//текстовый файл

function vvod_i(i,j:integer): Real;
begin
Q:= Fi/Ht;

if (i <= 1) then result := 0 else

if (j <= 1) then result := 0 else

if (j = m) then result := 0 else

if (i <= Q) and (j = Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + ( (f *(Ht * Ht)) + (2 * Arr[i,j] - Arr[i-1,j])) else

if (i <= Q) and (j > Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i = Q) and (j < Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i >= Q) and (j > Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i >= Q) and (j < Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i >= Q) and (j = Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else;

end;

begin
write('Введите кол-во узлов n');Readln(n);
write('Введите кол-во узлов m');Readln(m);
Write('Введите плотность струны');Readln(a);
Write('Введите силу воздействующую на узел (в ньютонах)');Readln(F);
Write('Введите шаг по времени');Readln(Ht);
Write('Введите шаг между узлами');Readln(Hx);
Write('Введите точку приложения силы');Readln(ft);
write ('введите время действия силы'); Readln(fi);
//Заполняем матрицу очевидными значениями

Arr[0, 0] := 0;
for i := 0 to n do
begin
Arr[i, 0] := i;
end;
for j := 1 to m do
begin
Arr[0, j] := j;
end;
//теперь расчитывем значения в узлах от 1 до размерности
for i := 1 to n do

for j := 1 to m do

Arr[i, j] := vvod_i(i,j);
//посмотрим что получилось
assign(fCSV,'c:\temp\out1.txt');//связываем файловую переменную с файлом
rewrite(fCSV);//открываем файл на перезапись
for i := 0 to n do
begin
for j := 0 to m do
begin
//выводим на экран
write(Arr[i, j]);
write(chr(09));
//выводим в файл
write(fCSV,Arr[i, j]);
write(fCSV,';');
end;
writeln;
//переход на новую сртоку
write(fCSV,chr(13));
write(fCSV,chr(10));
end;
close(fCSV);//закрываем файл
assign(s,'c:\temp\out2.txt');//связываем файловую переменную с файлом
rewrite(s);//открываем файл на перезапись
WriteLn(s, 'Plot ''do'' with it pause 0,5');
Close(s);
end.

Меня очень интересует нормально-ли идет расчет Q, и не изменится-ли рез-тат если меня местами уравнение(вместе с ограничениями) т.е. зависит-ли от места, то как будет идти расчет?

Как видите, расставлено много ограничений, но похоже он не верно видит массив...

Ошибка похоже где-то здесь (ниже),по тому как до сего момента всё в порядке

if (i <= Q) and (j > Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i = Q) and (j < Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i >= Q) and (j > Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i >= Q) and (j < Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else

if (i >= Q) and (j = Ft) then

result := (((a * a) * (Arr[i,j+1] - (2 * Arr[i,j]) + Arr[i,j-1]) * (Ht*Ht))/(Hx*Hx)) + (2* Arr[i,j]) - Arr[i-1,j] else;