Моделирование систем

 

Этап 1. Содержательная постановка задачи

 

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

 

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

 

 Рассмотрим одну такую простую  модель. Эта модель будет полезна  для прогноза событий и связанных с ними социально-экономических процессов.

 

Этап 2. Формулировка гипотез, построение, исследование модели

 

 Структура производства и  сбыта часто зависит от изменений  в окружающей среде (социально-экономических  условий).

 

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

 

 Пусть x(t) - величина основных  фондов (капитала) в момент времени  t, где 0 £ t £ N. Через промежуток  времени Dt она будет равна x(t+Dt).

 

Абсолютный прирост равен 

 

 

Относительный прирост будет равен 

 

 

Примем следующие гипотезы:

социально-экономические условия  производства достаточно хорошие и  способствуют росту производства, а  поток инвестиций задается в виде известной функции y(t);

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

 

 

 Считая Dt ® 0, а также учитывая  определение производной, получим  из предыдущего соотношения следующее  математическое выражение закона  изменения величины капитала - математическую модель (уравнение) динамики капитала (такие уравнения называются дифференциальными):

 

 

где х(0) - начальное значение капитала в момент времени t=0.

 

 Эта величина х0 везде в  дальнейшем будет считаться заданной. Эта простейшая модель динамики величины капитала.

 

 Эта простейшая модель не  отражает того факта, что социально-  экономические ресурсы производства  таковы, что между выделением  инвестиций и их введением  и использованием в выпуске  новой продукции проходит некоторое  время - лаг. Учитывая это можно записать модель (1) в виде:

 

 

 Данной непрерывной, дифференциальной, динамической модели можно поставить  в соответствие простую дискретную  модель:

 

 

где n - предельное значение момента  времени при моделировании. Эта  дискретная модель получается из непрерывной при Dt=1, а также заменой производной x'(t) на относительное приращение Dt (замена, как это следует из определения производной, справедлива при малых Dt).

 

Этап 3. Построение алгоритма и программы  моделирования 

 

 Рассмотрим  для простоты режим моделирования когда m, c, y - известны и постоянны, а также рассмотрим наиболее простой алгоритм моделирования в укрупнённых шагах.

Ввод входных  данных для моделирования:

с=х(0) - начальный  капитал;

n - конечное  время моделирования; 

m - коэффициент амортизации;

s - единица  измерения времени; 

y - инвестиции.

Вычисление xi от i=1 до i=n по рекуррентной формуле (2).

Поиск стационарного  состояния - такого момента времени j, 0<=j<=n, начиная с которого все  хj, хj+1,...,хn постоянны или изменяются на малую допустимую величину e>0.

Выдача  результатов моделирования и, по желанию пользователя, графика.

 

 Приведём  программу на Паскале для имитационного  моделирования (компактность программы  «принесена в жертву» структурированности).

 

 

 

PROGRAM MODFOND;

{Исходные  данные находятся в файле in.dat текущего каталога}

{Результаты  записываются в файл out.dat текущего  каталога}

Uses

  Crt, Graph, Textwin;

Type

  Vector = Array[0..2000] of Real;

  Mas    = Array[0..2000] of LongInt;

Var

  Time,Lag,t,dv,mv,i,yi,p           : Integer;

  tmax,tmin                         : LongInt;

  a,b,m,X0,maxx,minx,aa,bb,cc,sx,tk : Real;

  x                                 : Vector;

  ax,ay                             : Mas;

  ch                                : Char;

  f1,f2                             : Text;

{----------------------------------------------------------------------------}

Procedure InputKeyboard; { Ввод с клавиатуры }

Begin

OpenWindow(10,5,70,20,' Ввод  данных ',14,4);

 ClrScr;  WriteLn;  WriteLn('Введите время Т прогнозирования системы:');

Repeat

   Writeln('Для  удобства построения графика  введите Т не меньше 2');

   Write('Т =');     ReadLn(Time);

until Time>=2;

WriteLn('Введите лаг:');

 Repeat

   Write('Лаг  должен быть строго меньше  Т  - ');   ReadLn(Lag);

 until Lag0 )  a = ');  ReadLn(a);

Write('( b>0 )  b = ');  ReadLn(b);

Repeat

   Write('Введите коэффициент амортизации ( 0 < M < 1 ) - ');  Readln(m);

 until (m<1) and (m>0);

Write('Введите значение фондов  в начальный момент - ');  Readln(X0);

 CloseWindow;

end;

{----------------------------------------------------------------------------}

Procedure InputFile; { Ввод из файла }

Begin

Assign(f1,'in.dat'); Reset(f1);

Readln(f1,time,lag,a,b,m,X0); Close(f1);

End;

{----------------------------------------------------------------------------}

Procedure OutputFile; { Запись результата  работы в файл }

Begin

Assign(f2,'out.dat');  Rewrite(f2);

 WriteLn(f2,'   Результаты моделирования:');

WriteLn(f2,'Значение фондов в заданное  время Т = ',x[time]:4:2);

WriteLn(f2,'Максимальное значение фондов = ',maxx:4:2);

Write(f2,'Минимальное значение фондов = ',minx:4:2);

Close(f2);

End;

{----------------------------------------------------------------------------}

Procedure InputRnd; { Ввод случайными числами  }

Begin

Randomize;

Repeat   Time:=Random(90);  until Time>=2;

Repeat    Lag:=Random(80);   until Lagmaxx

    then  begin  maxx:=x[t];  tmax:=t;  end

    else if x[t]');  OutTextXY(67,47,'^');  OutTextXY(57,415,'0');

 OutTextXY(80,45,'X(T) - (Величина основных фондов производства)');

 OutTextXY(590,415,'T');  OutTextXY(540,430,'(Время)');  SetColor(2);

For i:=1 to 16 do Line(67,70+20*i,70,70+20*i);

For i:=1 to 24 do Line(70+20*i,410,70+20*i,413);

Mas_OX; Mas_OY;

For t:=0 to time do

  Вegin

   SetColor(Blue);  Circle(ax[t],ay[t],2);

   SetFillStyle(SolidFill,Red); FloodFill(ax[t],ay[t],Blue);

  End;

SetColor(Red);  SetLineStyle(3,1,1);

Line(70,ay[time],ax[time],ay[time]);

Line(ax[time],ay[time],ax[time],410);

Ipol(0,x[0],1,x[1],2,x[2]);

For i:=ax[0] to ax[2] do

  Begin

   sx:=p*(i-70)/20;

   yi:=410-round(20*(aa*sx*sx+bb*sx+cc)/tk);

   SetColor(Red);  Circle(i,yi,1);

  End;

For t:=1 to Time-2 do

  Begin

   Ipol(t,x[t],t+1,x[t+1],t+2,x[t+2]);

   For i:=ax[t+1] to ax[t+2] do

    Begin

     sx:=p*(i-70)/20;

     yi:=410-round(20*(aa*sx*sx+bb*sx+cc)/tk);

     SetColor(Red); Circle(i,yi,1);

    End;

  End;

ReadKey; CloseGraph;

End;

{----------------------------------------------------------------------------}

Begin

While true do

  Begin

   ClrScr;  TextBackGround(2);  Window(1,1,80,25);  ClrScr;

   OpenWindow(30,22,50,24,' Нажмите клавишу... ',4,1);

   OpenWindow(5,5,75,16,' Динамика фондов производства ',14,5);

   ClrScr;  WriteLn;

   WriteLn('  Пусть х(t) - основные  фонды в момент времени t,  y(t) -');

   WriteLn(' инвестиции,  m - коэффициент  амортизации фондов.');

   WriteLn(' Модель динамики основных  фондов (L - лаг):');

   Write('    x`(t) = y(t-L) - mx(t),  где х(0) = Хо,  y(t)=at+b, ( a,b>0 ).');

   ReadKey;  CloseWindow;

   OpenWindow(15,10,65,17,' Выбирите  вариант  ввода: ',15,0);

   ClrScr;   WriteLn;

   WriteLn('        С  клавиатуры        - < 1 >');

   WriteLn('        Из файла            - < 2 >');

   WriteLn('        Случайными  числами  - < 3 >');

   WriteLn('        Выход               - ');

   ch:=ReadKey;

   Сase ch of

     #49:  InputKeyboard;

     #50:  Вegin InputFile; OutputScreen; Еnd;

     #51:  Вegin InputRnd; OutputScreen; End;

     #27:  Halt(1);

    End;

   CloseWindow; Worker; OutputFile; OpenWindow(22,10,58,14,'',15,5);

   ClrScr;   WriteLn;

   Write('Для просмотра графика нажмите ввод');  ch:=ReadKey;

   If ch=#13 then begin Graf; RestoreCrtMode; end;

   CloseWindow;  TextBackGround(15);  Window(1,1,80,25);

   ClrScr; OpenWindow(15,10,65,16,'',15,6); ClrScr;  WriteLn;

   WriteLn('          Хотите еще моделировать ?');   WriteLn;

   WriteLn('Для выхода нажмите          -        < Esc >');

   WriteLn('Для продолжения нажмите любую другую клавишу');

   ch:=ReadKey;

   If ch=#27  then Halt(1);

   CloseWindow;

  End;

ClrScr;  TextBackGround(0);

End.

 

Этап 4. Проведение вычислительных экспериментов 

 

 Эксперимент 1. Поток инвестиций  постоянный и в каждый момент  времени равен 111. В начальный момент капитал - 1000 руб. Коэффициент амортизации - 0.0025. Построить модель динамики (посуточно) и найти величину основных фондов через 50 суток, если лаг равен 10 суток.

 

 Эксперимент 2. Основные фонды  в момент времени t=0 была равны  50000. Через какое время общая их сумма превысит 1200000 руб., если поток инвестиций постоянный и равен 200, а m=0.02, T=5?

 

Этап 5. Модификация (развитие) модели

 

 Для модели динамики фондов  с переменным законом потока  инвестиций:

а) построить гипотезы, модель и алгоритм моделирования;

б) сформулировать планы вычислительных, компьютерных экспериментов по модели;

в) реализовать алгоритм и планы  экспериментов на ЭВМ.

 

Компьютерное моделирование и  вычислительные эксперименты

 

 Рассмотрим следующую модель течения эпидемии. Пусть существует группа из n контактирующих индивидуумов, в которой в момент времени t имеется x восприимчивых индивидуумов, y - источников инфекции и z - удаленных (изолированных или выздоровевших и ставших невосприимчивыми к инфекции) индивидуумов. Таким образом, имеем: x+y+z=n.

 

 Определим исходные параметры  системы: b — частота контактов  между членами группы; g — частота  случаев удаления; m — скорость  пополнения восприимчивых индивидуумов  извне (она же и есть скорость  гибели индивидуумов, удаленных из популяции); y* — критическое значение, при котором начинается эпидемия.

 

 После какой — либо одной  вспышки эпидемии, в результате  которой плотность восприимчивых  индивидуумов упадет ниже критического  значения y*, наступает период относительного затишья, длящийся до тех пор, пока снова не будет достигнуто критическое значение y* и не возникнет новая вспышка. За время Dt группа восприимчивых индивидуумов, с одной стороны, уменьшается на b xy Dt за счет заражения части из них, а с другой — увеличивается на n Dt.

 

 При t=0 заданы x(0)=x0, y(0)=y0, z(0)=z0 - число  восприимчивых, инфекционных и  удаленных индивидуумов соответственно.

 

 Рассматривается промежуток  времени t (0 £ t £ T). При малых  Dt=1 (например, минута) получим соотношения  (дискретную модель):

 

 

Приведём результаты проведённых  вычислительных экспериментов.

 

 

 

 

 

 

 Первые два случая демонстрируют  периодичность эпидемии. В последнем  примере имеет место случай  стационарного состояния, т.е.  после некоторого момента времени t количества индивидуумов в группах x, y, z остаются постоянными.

 

 Программу компьютерного моделирования  и вычислительные эксперименты  реализовали студентки математического  факультета КБГУ Буздова Асият  и Мильман Евгения. 

 

Компьютерное моделирование.

 

 Рассмотрим проблему расчета  влажности почвы с учетом накапливаемой  биомассы и прогнозирования урожайности  сельхозкультур по заданной (экологически  обоснованной) влагообеспеченности  корнеобитаемого слоя почвы. Разработать  соответствующую компьютерную моделирующую среду, которая позволяет решать задачи прогноза влажности корнеобитаемого слоя почвы и урожайности (биомассы) сельхозкультур на заданный момент времени с развитыми интерфейсными средствами, рассчитанными на неподготовленного пользователя-агронома, эколога и др. Опишем одну такую программную среду, реализованную реально в среде Delphi 2.0 Windows 95 автором и студентами КБГУ Кирьязевой С.К. и Кирьязевым Д.А. для расчета влажности почвы и определения урожайности сельхозкультур. Общение с пользователем осуществляется с помощью диалогового окна “Расчет влажности и урожайности”, содержащего 5 страниц: “Эксперимент”, “С/х культура”, “Регион”, “Рабочая” и “Результат”.

 

 Страница “С/х культура”-  для ввода входной информации  по культуре.

 

 

Страница “Регион” - для ввода информации по региону эксперимента.

 

 

Страница “Эксперимент” выглядит следующим образом.

 

 

 Данная страница предназначена  для ввода названия эксперимента, выбора названия культуры из  списка имеющихся (которые были  введены на странице “С/х культура”), региона проведения эксперимента из списка регионов, введенных на странице “Регион”, ввода даты посева культуры и даты снятия урожая, для ввода данных по величинам, зависящим от вида культуры, типа почв (по фазам вегетации). Предусмотрена возможность выбора: вычислять ли значения коэффициентов e и l или же рассматривать их как постоянные величины? После заполнения страницы “Эксперимент”, можно произвести расчет влажности почвы и прогноз урожайности культуры. Для этого достаточно лишь нажать кнопку “Произвести расчет” на странице “Эксперимент”. После этого автоматически раскрывается страница “Результат” с таблицей рассчитанных величин и выводится график зависимости влажности почвы от времени (1 - синим цветом), накопления биомассы растений от времени при вычисленной влажности (2 - зеленым цветом) и оптимального развития растений - по экспериментальным данным за прошлый год (3 - красным цветом).