Программирование математических приложений

В данной статье рассматриваются основные аспекты программирования пользовательских математических приложений. Будут описаны общие принципы и приближенный алгоритм действий по распознаванию формулы, построению 3D графиков и нахождению производной формулы. Прежде, чем вы начнете читать статью, посмотрите программу Математика v0.5 (автор программы Дьяченко Евгений), чтобы иметь представление как все должно выглядеть. Статья состоит из трех самостоятельных частей, которые можно читать в отдельности друг от друга.

1. Распознавание формулы из строки

Главная проблема, которая встает перед программистом, это написать программу так, чтобы для каждой новой формулы не приходилось ее перекомпилировать. Программа должна считывать формулу из какой-либо строки (например, form1.edit1.text) и затем ее обрабатывать. Для того, что бы рассчитывать строку необходимо также иметь массив значений переменных величин, которые будут использоваться в процессе распознавания. При записи формулы пользователь должен иметь возможность использовать стандартные математические функции Delphi, с их стандартными обозначениями. Формула может состоять из любых комбинаций функций и арифметических действий, включая вложенные функции, если это позволяют законы математики.

Итак, опишем процесс распознавания. Вначале, в целях повышения быстродействия программы, необходимо укоротить строку. Для этого заменим все стандартные обозначения математических функций, другими обозначениями, состоящими из одного символа (например, cos на С, sin на S). Важно проследить, чтобы обозначения не дублировались. Также можно вырезать из строки все пробелы.

Следующим шагом будет поиск в строке переменных и замена переменных их значениями, которые берутся из вышеописанного массива значений переменных. Для этого вырезаем из строки название переменной и вставляем значение переменной, которое с помощью функции floattostr превратили в строку. Функция floattostr превращает число в строку, при этом в получившейся строке, число имеет очень много знаков, что замедляет процесс распознавания, т.к. увеличивает длину исходной строки. Чтобы этого избежать, можно использовать вместо floattostr функцию floattostrf , которая имеет дополнительные параметры для управления процессом преобразования числа в строку. Подробнее о функции читайте в хелпе Delphi.

Теперь нам необходимо в строке с формулой найти самую глубокую вложенную функцию. Оказывается, это легко осуществить, если проссматривать строку с конца циклом

for i:=length(stroka) downto 1 do …

Первая попавшаяся функция и будет самой глубокой. Берем аргумент (он может содержать арифметические действия) этой функции и рассчитываем его специальной подпрограммой, которая рассчитывает формулы, не содержащие математических функций. Данная подпрограмма будет описана ниже. Вставляем вместо сложного аргумента нашей функции значение (число) этого аргумента. Сразу после этого заменяем функцию ее значением. Т.е. с помощью оператора

case of

сравниваем название нашей функции с названиями всех известных математических функций, и, если есть совпадение названий, с помощью стандартной подпрограммы strtofloat переводим аргумент функции в число, рассчитываем соответствующую функцию от этого аргумента и заменяем обозначение функции в исходной строке ее значением.

Далее продолжаем проссматривать строку все тем же циклом, при этом все функции будут автоматически правильно заменяться их значениями. В конце останется строка, не содержащая ни одной математической функции. Рассчитаем эту строку подпрограммой, которая упоминалась выше и будет описана в следующем абзаце.

Итак, нам необходима подпрограмма, которая заменяет строку, содержащую арифметические действия и приоритетные скобки, на результат расчета этой строки.

Вначале предположим, что в строке нет приоритетных скобок, только арифметические действия. Программа вначале ищет во всей строке арифметическое действие возведения в степень, когда находит, берет число слева и возводит в степень числа справа, затем продолжает поиск, пока не достигнет конца. Затем проделывает тоже самое сначала для деления и умножения, потом для сложения и вычитания. Этим мы обеспечим принятый в математике приоритет арифметических действий.

Теперь рассмотрим случай, когда в строке есть приоритетные скобки. Программа должна последовательно находить самые глубокие скобки и производить над ними то, что было описано абзацом выше. Этот алгоритм легче представить кодом:

s1:='';
for i:=length(a) downto 1 do //просматривает строку а с конца 
begin
 if a[i]=')' then ni:=i; //запоминает положение последней ) 
 if a[i]='(' then //если встречает (, значит это и есть самая глубокая ( ) 
 begin
 for j:=i+1 to ni-1 do s1:=s1+a[j];
 delete(a,i,2+length(s1)); //удаляет исходную скобку… 
 insert(rasstr(s1),a,i); //… и вставляет вместо нее значение скобки
 //возвращается выше на скобку и продолжает далее цикл 
 for j:=i to length(a) do 
 if a[j]=')' then begin ni:=j; break; end;
 s1:='';
 end;
end;
a:=rasstr(a); //расчет конечного результата 

Здесь rassrt - подпрограмма для решения строчек, в которых нет приоритетных скобок (только арифметические действия).

Последний раз rasstr применяется для решения конечной строки, из которой мы убрали все приоритетные скобки. Переводим строку в число, и получаем значение исходной строки с формулой, т.е. результат. Хочу предупредить о возможных проблемах:

1. При умножении чисел часто возникает ситуация, когда второй множитель отрицательный - программа не должна путать его со знаком вычитания. Необходимо переносить минус в начало члена, т.е. изменять знак перед членом.

2. Числа в компьютере могут представляться в экспонентальном виде (например, 1Е-5, что равно 0.00001). Может возникнуть проблема с минусом в показателе. Необходимо или особо описать эту ситуацию в программе, или использовать для преобразования чисел в строку функцию FloatToStrF (о функции смотрите в справочнике Delphi). Первый путь предпочтительнее, так как он будет быстрее работать (в большинстве случаев).

3. В Delphi нет функции возведения произвольного числа в произвольную степень. Это решается следующим образом: a^b=exp(ln(a)*b).

2. Построение трехмерных графиков

Необходимо отметить, что графики в приведенном алгоритме рисуются стандартными процедурами LineTo и MoveTo, т.е. данный алгоритм вы можете реализовать в программе, не используя стандартные процедуры Delphi по работе с трехмерным пространством.

Прежде всего, опишу, как можно представить трехмерные графики в памяти компьютера. Допустим, нам надо построить какой-либо трехмерный график на диапазоне по Ox и Oy от -10 до 10 с точностью в 100 шагов на этом диапазоне. Возьмем двухмерный массив

array [1..100,1..100] of real;

Обозначим первую координату массива за переменную х, второю за y, значение ячейки [x,y] будет значением переменной z. Разделим диапазон на 100 шагов, для каждого шага будем рассчитывать z=f(x,y), и заносить значение z в ячейку массива [<номер шага по х>,<номер шага по y>]. Т.е. двумя циклами (второй цикл вложенный) просматриваем все комбинации переменных x и y с определенным шагом, и для каждой комбинации заполняем соответствующую ячейку массива.

Для вращения графика используется изменение трех переменных: alf - поворот вокруг оси Oz, bet - поворот вокруг оси Oy, gam - поворот вокруг оси Ox. При вращении графика наш массив, содержащий график остается неизменным, мы просто под другими углами проецируем его на поверхность экрана.

Переменные alf, bet, gam изменяются при нажатии на поле рисования и передвижения мыши. Это происходит следующим образом:

procedure TForm1.Image1MouseDown(Sender: TObject; Button: TMouseButton;
 Shift: TShiftState; X, Y: Integer);
begin 
 md:=true; xa:=x; ya:=y; 
end;
procedure TForm1.Image1MouseUp(Sender: TObject; Button: TMouseButton;
 Shift: TShiftState; X, Y: Integer);
begin 
 md:=false; xa:=x; ya:=y; 
end;
procedure TForm1.Image1MouseMove(Sender: TObject; Shift: TShiftState; X,
 Y: Integer);
begin
 if md then 
 begin 
 alf:=alf-(x-xa)*2*pi/300; //число 300 означает, на сколько пикселей надо 
 bet:=bet-(y-ya)*2*pi/300; //передвинуть мышь для поворота графика на 360 градусов 
 xa:=x; ya:=y;
 risovanie; // процедура для перерисовки графика, с новыми значениями углов. 
 end;
end;

В процедуре, которую мы назвали risovanie, строятся линии специальными подпрограммами LineToXYX и MoveToXYZ, которые будут описаны ниже. Для этого двумя циклами (второй вложенный) просматриваются все ячейки массива и аргументам функций LineToXYX и MoveToXYZ присваиваются значения (<начало диапазона х>+<длина диапазона х>*<номер текущего шага по х>/<количество шагов по х>, <тоже самое по у>, <значение массива в точке [<шаг по х>, <шаг по у>]>). Сначала строится одна половина решетки, потом другая. Это делается двумя последовательными циклами.

Создадим две процедуры MoveToXYZ и LineToXYZ для проецирования линий из трехмерного пространства на поверхность экрана.

procedure LineToXYZ(x,y,z:real);
begin
 form1.image1.canvas.lineto(trunc(((y*cos(alf)-x*sin(alf))*cos(gam)+(z*cos(bet)
 +(x*cos(alf)+y*sin(alf))*sin(bet))*sin(gam))*zum+form1.image1.Width/2),
 trunc(-((z*cos(bet)+(x*cos(alf)+y*sin(alf))*sin(bet))*cos(gam)-(y*cos(alf)-
 x*sin(alf))*sin(gam))*zum+form1.Image1.Height/2));
end;
procedure MoveToXYZ(x,y,z:real);
begin
 form1.image1.canvas.moveto(trunc(((y*cos(alf)-x*sin(alf))*cos(gam)+(z*cos(bet)
 +(x*cos(alf)+y*sin(alf))*sin(bet))*sin(gam))*zum+form1.image1.Width/2),
 trunc (-((z*cos(bet)+(x*cos(alf)+y*sin(alf))*sin(bet))*cos(gam)-(y*cos(alf)-
 x*sin(alf))*sin(gam))*zum+form1.Image1.Height/2));
end;

Переменные alf, bet, gam - глобальные, zum - переменная, которая определяет коэффициент увеличения графика, тоже глобальная. Trunc - стандартная функция, возвращает значение целой части аргумента и имеет значение типа Integer. График рисуется относительно центра поля рисования.

Хочу отметить, что поворот вокруг оси Ox не добавляет наглядности графику. Если его исключить, формула сократиться вдвое, следовательно, повыситься скорость работы. Вот как выглядят сокращенные процедуры:

procedure LineToXYZ(x,y,z:real);
begin
 form1.image1.canvas.lineto(trunc((y*cos(alf)-x*sin(alf))*zum+form1.Image1.width/2),
 trunc(-(z*cos(bet)+(x*cos(alf)+y*sin(alf))*sin(bet))*zum+form1.image1.height/2));
end;
procedure MoveToXYZ(x,y,z:real);
begin
 form1.image1.canvas.moveto(trunc((y*cos(alf)-x*sin(alf))*zum+form1.Image1.width/2),
 trunc(-(z*cos(bet)+(x*cos(alf)+y*sin(alf))*sin(bet))*zum+form1.image1.height/2));
end;

3. Нахождение производной аналитическим методом

Эта задача не очень сложная (впрочем, как и все остальные), но достаточно трудоемкая. Приведу словесный алгоритм выполнения этой операции.

1. Замена стандартных обозначений функций их краткими служебными названиями для ускорения работы программы.
2. Повторять замену производных суммы (разности) и произведения (частного) по соответствующим формулам, пока не останется ни одной производной суммы (разности) и произведения (частного). При этом каждый раз необходимо раскрывать скобки.
3. Замена всех производных их табличными значениями.
4. Повторять пункты 2 и 3, пока не останется ни одной производной.
5. Упростить полученное выражение: выкинуть из произведения все единицы, выкинуть все члены, в которых встречаются нули и т.п.
После каждой замены необходимо производить раскрытие скобок, без этого программа не будет работать. Код, который демонстрирует словесный алгоритм:

function proizvod(a:string):string; //общая производная 
var a:string;
begin
 repeat
 repeat
 zamen:=false;
 a:=psum(a); //производная суммы и разности 
 a:=pproi(a); //производная поизведения и частного 
 until zamen=false;
 a:=ptab(a); //производная табличных величин 
 until zamen=false;
 a:=pupro(a); //упрощение внешнего вида формулы 
 proizvod:=a;
end;

Переменной zamen (типа Boolean) после любой замены, произошедшей в формуле, присваивается значение true.

Для замены производной суммы, разности, произведения и частного подпрограммы очень похожи, поэтому я объединю их описание. В строке производится поиск всех производных. Как только находиться какая-нибудь производная, она проверяется на наличие внутри символов арифметических действий. Если они есть, то производная раскладывается по соответствующей формуле, и переменной zamen, о которой говорилось выше, присваивается значение true. Далее продолжается поиск в строке других производных. По завершению поиска необходимо произвести раскрытие приоритетных скобок. Алгоритм будет быстрее работать, начать поиск в строке с конца. При каждой замене строка удлиняется, если производить поиск с конца, программа не будет заново просматривать то, что она только что вставила в строку.

Замена производных табличных величин происходит аналогично. Это сделать очень легко, т.к. известно, что в производных нет знаков арифметических действий.

Вообще-то, упрощение внешнего вида формулы можно не производить, но тогда формула будет не наглядной для человека, и будет медленно обрабатываться компьютером. Поскольку эта подпрограмма не является принципиально важной, опущу ее описание.

Заключение

В заключении, хочу отметить, что приведенный алгоритм является не совсем полным. Он будет работать с определенными ограничениями. Для того, что бы программа умела работать с абсолютно любыми формулами, необходимо решить огромное количество мелких и достаточно простых проблем, их можно назвать предельными случаями, которые требуют отдельного описания в программе. Основные проблемы возникают со знаком минуса и раскрытием скобок. В процессе написания программы эти проблемы быстро решаются, если их обнаружить, что не всегда легко. Поэтому, советую использовать систему автоматического тестирования программы на выполнение математических операций (функций, действий, историй). Тесты должны запускаться после каждого изменения в программе. В них необходимо описать общие случаи и предельные, для которых заранее известен результат.

Если у вас есть замечания по статье или предложения пишите на e-mail: ka82507@mail.tomsknet.ru. Связаться со мной также можно по тел. 382-2-54-13-60. Дьяченко Евгений.