Лабораторный практикум по моделированию

МИНОБРНАУКИ РОССИ Федеральное государственное бюджетное образовательное  учреждение высшего профессионального образования «Чувашский государственный университет имени И.Н. Ульянова» В. К. НИКИШЕВ Лабораторный практикум  по моделированию  Учебное  пособие для студентов кафедры компьютерных технологий Чебоксары, 2013 г Печатается по решению ученого совета ГОУ ВПО «Чувашский   государственный   университет  им. И. Н. Ульянова»  Рецензенты: Ванюлин   А.   Н.   ­  канд.   тех.   наук,   доцент кафедры   информационных   технологий   Чебоксарского кооперативного   института   Московского   университета потребительской кооперации.                    Зайцева В. П. ­   канд. фил. наук, доцент кафедры Чебоксарского информационных педагогического университета им. И. Я. Яковлева     технологий .               В пособии представлены примеры выполнения лабораторных работ по моделированию на языках программирования   Delphi, VS 2010 ( VC#,VC++, VB.NET) и в среде моделирования  SciLab (MatLab). Тематика лабораторных работ соответствует рабочей программы   по   моделированию,   которая   составлена   в соответствии   со   стандартом   образования.   Каждая   работа содержит     условие   задачи,   алгоритм   в     виде   блок­схемы, программы на одном из языков программирования  и результаты исследований.   Задания   для   выполнения   лабораторных   работ представлены   в   конце   учебного   пособия.   Лабораторный практикум   могут   использовать   преподаватели,   магистры   и аспиранты,   которые   изучают   или   преподают   основы математического   моделирования,   особенно   с   реализацией   на языках программирования VC#, VC++                  ©  Никишев В. К., 2013 2 Оглавление 1.  Введение 2. Лабораторные работы по моделированию…………………4  2.1  Лабораторная работа №1      Тема: Методы исследования дифференциальных     уравнений  с использованием разных информационных технологий…… 8  2.2  Лабораторная работа № 2  Тема: Свободное падение тел…………………………………   8  2.3  Лабораторная работа № 3  Тема: Исследование динамики объектов,  брошенных  под углом к горизонту………………………………………….40     2.4  Лабораторная работа №4                              Тема: Исследования динамики  популяций хищника и жертвы………………………………….68    2.5 Лабораторная работа №5                                   Тема: Моделирование движения  заряженных частиц и небесных тел……………………………81     2.6   Лабораторная работа № 6                            Тема: Исследование динамики полета ракеты.……..………..98     2.7  Лабораторная работа № 6                         Тема: Разработка информационной  модели студента……..115 2.8  Лабораторная работа № 8                      Тема: Исследование экономических систем ……..………….127 Тема: Исследование  систем массового обслуживания …     143  2.9  Лабораторная работа № 9                      Тема: Исследование оптимальных систем……………………150  2.10 Лабораторная работа №10         Тема.  Моделирование объектов методом                        пространства состояния…………………………………   159    Литература   3 Введение. В   настоящее   время   большое   внимание   уделяется   вопросам моделирования различных систем с использованием современных языков   программирования   (   Visual   Basic,   Delphi,    VC#,VC++, VB.NET  ) и информационных программ, как, например, Excel, MathCad, MatLab, Maple,  SciLab.  В отличии от программирования, где разрабатываются алгоритм и   программа   для   решения   какой   либо   задачи   для   получения результата   решения   при   заданных   исходных   данных,   в моделировании   разрабатываются   алгоритм   и   программа   для исследования   систем,   объектов   или   процессов.   Необходимо помнить,   что   моделирование   ­   это   исследование   систем,   это вычислительный   эксперимент.   А   исследование   обычно проводится с учетом воздействия на модель, представленной в математической   или   иной   формах,   различных   входных параметров или изменение различных коэффициентов, которые входят   в   уравнение   модели.   В   результате   проведения вычислительного   эксперимента   по   полученным   результатам можно   сделать   соответствующие   выводы   по   устойчивости систем, точности систем, управлению объектов или в целом по работе  какой ­  либо  информационной  системы  при различных воздействиях   на   систему.   Поэтому   в   отличие   от   простой программы   необходимо   разработать   проект   для   исследования системы.   Такой   проект   может   иметь   следующую   структуру: получение   результата   моделирования   при   конкретных параметрах,   при   изменении   параметров   в   определенных интервалах и получения так называемого среза результата при изменении исследуемого параметра в определенных интервалах. 4 2. Лабораторные работы по моделированию   2.1  Лабораторная работа №1      Тема:   Методы   исследования   объектов,   динамика   которых описывается дифференциальными   уравнениями с   использованием   языков   программирования  Delphi,  VC+ +.NET, VC#.  NET    Цель занятия: 1   Получить   практические   навыки   исследования   систем (объектов), динамика которых описывается дифференциальными уравнениями 1­го порядка  2.   Научиться   разрабатывать   алгоритм   и   программу   с использованием   языков   программирования  Delphi,  VC++.NET, VC#.NET  3. Практически усвоить численные методы Эйлера и Рунге­Кутта для  решения дифференциальных уравнений 1 –го порядка. Задачи  занятия: 1. Разработка алгоритма в виде блок ­ схемы 2. Построение графиков  кривых y(x), dy/dx при параметрах a­ const и var,  3. Анализ результатов исследований. Модели объектов исследования Программа исследования aӳ+bý + cy=f 5 1. a,c,f ­ const, t – var ( t0 – tk, h= 0.1, 0.01) 2. a,f ­ const, t ­var( t0 – tk, h= 0.1, 0.01), b­ var Пример 1 Условие задачи: Составить алгоритм и проект моделирования  объекта, динамика которого описывается   дифференциальным    уравнением 1 –го порядка методом Эйлера                                    5 ý + 3 y=4                              Блок­схема алгоритма 6 Отчет по лабораторной работе 7 Титульный лист Форма для  исследования объекта Листинг программы ( язык программирования      Delphi) procedure TForm3.Button2Click(Sender: TObject); var: x:real; i:integer; function F(x,y:real):real; begin f:=5*x+(3*y)­4; end; begin 8 a:=StrToFloat(Edit1.text);  b:=StrToFloat(Edit2.text);  h:=StrToFloat(Edit3.text);  y1:=strtofloat(edit4.text);  //шаг memo1.lines.add(floattostr(x)); memo2.lines.add(floattostr(y)); x:=a;   Chart1.Series[0].clear; Chart1.Series[1].clear;  repeat   begin   y:=y1+h*f(x,y);    Chart1.Series[0].Add(y,FloatToStr(x),clRed);     Chart1.Series[1].Add(f(x,y),FloatToStr(x),clblue);    x:=x+h;   y1:=y;     memo1.lines.add(floattostrf(x,fffixed,7,3));     memo2.lines.add(floattostrf(y,fffixed,7,3));    end;   until x>b; end; Листинг программы ( язык программирования      VC++) Исходное дифференциальное уравнение ý ­ y ­ 2*sin(x) = 0; 9 private: System::Void button1_Click(System::Object^  sender,  System::EventArgs^  e) {  this­>Hide();  } private: System::Void button2_Click(System::Object^  sender,  System::EventArgs^  e) {  // int i = 1;  // int n = 100;   double f1,x, y1;   double h = Convert::ToDouble(textBox1­>Text);   double b = Convert::ToDouble(textBox2­>Text);   double y = Convert::ToDouble(textBox3­>Text);      x = 0;   y1 = y;   chart1­>Series[0]­>Points­>Clear();   chart1­>Series[1]­>Points­>Clear();   do   {   f1 = y ­ 2*sin(x);   x = x+h;   y = y+f1*h;    chart1­>Series[0]­>Points­>AddXY(x, y);   chart1­>Series[1]­ >Points­>AddXY(x,f1);    }   while (x<=b);      }};}};} 10 Язык программирования  C# private void button1_Click(object sender, EventArgs e)         {             double h = 0.01;             double y = 0;             double b = 2;             double x = 0;             double y1 = 4;             do             {                 double f = y + 2*Math.Sin(x);                 x = x + h;                 y = y1 + f * h;                 chart1.Series[0].Points.AddXY(x, y); 11 chart1.Series[1].Points.AddXY(x, f);                 y1 = y;             }             while (x <= b);         }                        // foreach (DataPoint p in chart1.Series[0].Points)            // {                //  Вывожу Х в лог                // textBox1.AppendText("X=" + p.XValue.ToString());                // textBox1.AppendText(Environment.NewLine);               //  listBox1.Items.Add("X=" + p.XValue.ToString());            // }             //   Y является массивом, поэтому пробегаю по массиву            // foreach (DataPoint yp in chart1.Series[0].Points)            // {                 //Вывожу Y                // textBox1.AppendText("Y=" + yp.ToString());                //  textBox1.AppendText(Environment.NewLine);                 // listBox1.AppendText("Y=" + yp.ToString());                 // listBox1.AppendText(Environment.NewLine);                // listBox2.Items.Add("Y=" + yp.ToString());            // }        // }               private void chart1_Click(object sender, EventArgs e)           {           }           private void button2_Click_1(object sender, EventArgs e)           {               Close();             }           }}   12 Пример 2.   Исследовать объект, динамика которого  представляется дифференциальным уравнением  ý ­ y ­ a*sin(x) = 0; (объектов),      Цель занятия: Получить практические навыки исследования систем     динамика   которых   описывается дифференциальными уравнениями  1­ го порядка при изменении параметра а  Задачи  занятия: 1. Разработка алгоритма в виде блок – схемы 2. Построение графиков  кривых y(x), dy/dx при параметрах a­var 13 private: System::Void button1_Click(System::Object^  sender,  System::EventArgs^  e) {  this­>Hide();  } private: System::Void button2_Click(System::Object^  sender,  System::EventArgs^  e) {  // int i = 1;  // int n = 100;   double f1,f2,x, y1;   double h = Convert::ToDouble(textBox1­>Text);   double b = Convert::ToDouble(textBox2­>Text);   double y2 = Convert::ToDouble(textBox3­>Text); 14 double a = Convert::ToDouble(textBox4­>Text);  double a1 = Convert::ToDouble(textBox5­>Text);      x = 0;  // double y1 = y2;   int a3=a;   //for int i=0;i<=1;i++)   //{   chart1­>Series[0]­>Points­>Clear();   chart1­>Series[1]­>Points­>Clear();    chart1­>Series[2]­>Points­>Clear();   chart2­>Series[0]­>Points­>Clear();    chart2­>Series[1]­>Points­>Clear();  chart2­>Series[2]­>Points­>Clear();   //}   do   {   // for (int i=1;i<=5;i++) //  {   a3=a;   f1 = y2 ­ a3*sin(x);  // int i=0; double  y = y2+f1*h; //int i=0;   chart1­>Series[0]­>Points­>AddXY(x, y);   chart2­>Series[0]­>Points­>AddXY(x,f1); a3=a3+a1; f1 = y2 ­ a3*sin(x);   y = y2+f1*h;  // f2 = y2 ­ 3*sin(x); //  y = y2+f2*h; // int i=1;   chart1­>Series[1]­>Points­>AddXY(x,y);      chart2­>Series[1]­>Points­>AddXY(x,f1);   a3=a3+a1; f1 = y2 ­ a3*sin(x);   y = y2+f1*h; 15 chart1­>Series[2]­>Points­>AddXY(x,y);    //y=y2;   chart2­>Series[2]­>Points­>AddXY(x,f1); x = x+h;   }   // a=a+3;   while (x<=b);       x=0;      } private: System::Void chart2_Click(System::Object^  sender,  System::EventArgs^  e) {  } private: System::Void textBox1_TextChanged(System::Object^   sender, System::EventArgs^  e) {  } }; } Лабораторная работа №2    Тема: Методы исследования объектов, динамика которых описывается дифференциальными  уравнениями 2­го порядка.    Цель занятия: Получить практические навыки исследования  систем (объектов), динамика которых описывается  дифференциальными уравнениями  2­ го порядка.  Задачи  занятия: 1. Разработка алгоритма в виде блок ­ схемы 2. Построение графиков  кривых y(x), dy/dx при параметрах a­ const,  3. Анализ результатов исследований. 16 Пример исследования на языке программирования VC++ 17 private: System::Void button1_Click(System::Object^  sender,  System::EventArgs^  e) {  this­>Hide();  } private: System::Void button2_Click(System::Object^  sender,  System::EventArgs^  e) {  // int i = 1;  // int n = 100;   double f1,x,y11,y22;   double h = Convert::ToDouble(textBox1­>Text);   double b = Convert::ToDouble(textBox2­>Text);   double y1 = Convert::ToDouble(textBox3­>Text);     double y2 = Convert::ToDouble(textBox4­>Text);      x = 0;   //y1 =0; //  y2=0;   chart1­>Series[0]­>Points­>Clear();   chart1­>Series[1]­>Points­>Clear();   do   { //  f1 = y ­ 2*sin(x);  x = x+h;  // y = y+f1*h;  y11=y1+h*y2;  y22=y2+h*(5­3*y2­4*y1);   //y = y1+h*(f1+(y ­ 2*x/y))/2;   chart1­>Series[0]­>Points­>AddXY(x, y11);   chart1­>Series[1]­>Points­>AddXY(x,y22);   y1=y11; 18 y2=y22;    }   while (x<=b);      } private: System::Void Form5_Load(System::Object^  sender,  System::EventArgs^  e) {  }};} Пример исследования на языке Делфи. Условие задачи: Составить алгоритм и проект моделирования  объекта, динамика которого описывается дифференциальным    уравнением 2 порядка методом Рунге Кутта                                                                         y’’+a*x*y’­y=0,4.                    Блок схема алгоритма 19 Программа на Delphi unit Unit2; interface var   Form2: TForm2; implementation uses Unit1,Unit4; procedure TForm2.SpeedButton1Click(Sender: TObject); begin form4.Show;form2.Hide; 20 end; procedure TForm2.SpeedButton2Click(Sender: TObject); var   i,p :integer; // переменные, использующиеся в циклах   a : integer; // а ­ параметр а в нашем уравнении   h,t:real;   // h ­ точность, t ­ шаг   n :integer;  //  число шагов  // x0,y0 :real; //  начальные значения x и y   tn, tk : integer; // границы   x,y,kx1,kx2,kx3,kx4,ky1,ky2,ky3,ky4,dx,dy :array[1..10000] of  real;   // x,y ­ значения x и y    // kx1­4,ky1­4 ­ переменные , принимающие участие в формуле  Рунге­Кутта   //dx, dy ­ отрезки x и y, из ,которых складывается решение  уравнения по Р­Кутту begin  a:=strtoint(edit1.Text);tn:=strtoint(edit2.Text); tk:=strtoint(edit3.Text);h:=strtofloat(edit4.Text); x[1]:=strtofloat(edit5.Text);y[1]:=strtofloat(edit6.Text); n:=trunc((tk­tn)/h); t:=tn; chart1.SeriesList[0].Clear;chart2.SeriesList[0].Clear; stringgrid1.RowCount:=n;   for i:= 1 to (n) do    begin  kx1[i]:=h*(­a*y[i]­6*x[i]);  kx2[i]:=h*(­a*y[i]­6*(x[i]+kx1[i]/2));  kx3[i]:=h*(­a*y[i]­6*(x[i]+kx2[i]/2)); kx4[i]:=h*(­a*y[i]­6*(x[i]+kx3[i]));  dx[i]:=1/6*(kx1[i]+2*kx2[i]+2*kx3[i]+kx4[i]);  x[i+1]:=x[i]+dx[i];  ky1[i]:=h*x[i];  ky2[i]:=h*(x[i]+ky1[i]/2);  ky3[i]:=h*(x[i]+ky2[i]/2);  ky4[i]:=h*(x[i]+ky3[i]); 21 dy[i]:=1/6*(ky1[i]+2*ky2[i]+2*ky3[i]+ky4[i]);  y[i+1]:=y[i]+dy[i];  t:=t+h;  stringgrid1.Cells[0,i]:=floattostr(t);  stringgrid1.Cells[1,i]:=floattostr(y[i]);  stringgrid1.Cells[2,i]:=floattostr(x[i]);  chart1.SeriesList[0].AddXY(t,y[i],'',clblue);  chart2.SeriesList[0].AddXY(t,x[i],'',clblue);    end;end; procedure TForm2.FormCreate(Sender: TObject); begin stringgrid1.Cells[0,0]:='x'; stringgrid1.Cells[1,0]:='f(x)'; stringgrid1.Cells[2,0]:='f`(x)';end;end. unit Unit5; interface     var   Form5: TForm5; implementation uses Unit2, Unit4; procedure TForm5.SpeedButton2Click(Sender: TObject); var   i,p :integer; // переменные, использующиеся в циклах   a,a1,da : real; // а ­ параметр а в нашем уравнении   h,t:real;   // h ­ точность, t ­ шаг   n,j :integer;  //  число шагов  // x0,y0 :real; //  начальные значения x и y   tn, tk : integer; // границы   x,y,kx1,kx2,kx3,kx4,ky1,ky2,ky3,ky4,dx,dy :array[1..10000] of  real;   // x,y ­ значения x и y на соотвествующем шаге   // kx1­4,ky1­4 ­ переменные опринимающие участие в формуле  Рунге­Кутта   //dx, dy ­ отрезки x и y из,которых складывается решение ур­ия  по Р­Кутту begin 22 for j:= 0 to 5 do begin chart1.SeriesList[j].Clear;chart2.SeriesList[j].Clear; end; a:=strtofloat(edit1.Text);a1:=strtofloat(edit7.Text); da:=strtofloat(edit8.Text); j:=0;a:=a­da; repeat tn:=strtoint(edit2.Text);tk:=strtoint(edit3.Text); h:=strtofloat(edit4.Text);x[1]:=strtofloat(edit5.Text); y[1]:=strtofloat(edit6.Text);n:=trunc((tk­tn)/h); t:=tn;chart1.SeriesList[0].Clear;chart2.SeriesList[0].Clear;   for i:= 1 to (n) do    begin  kx1[i]:=h*(­a*y[i]­6*x[i]);  kx2[i]:=h*(­a*y[i]­6*(x[i]+kx1[i]/2));  kx3[i]:=h*(­a*y[i]­6*(x[i]+kx2[i]/2));  kx4[i]:=h*(­a*y[i]­6*(x[i]+kx3[i]));  dx[i]:=1/6*(kx1[i]+2*kx2[i]+2*kx3[i]+kx4[i]);  x[i+1]:=x[i]+dx[i]; ky1[i]:=h*x[i];  ky2[i]:=h*(x[i]+ky1[i]/2); ky3[i]:=h*(x[i]+ky2[i]/2);  ky4[i]:=h*(x[i]+ky3[i]);  dy[i]:=1/6*(ky1[i]+2*ky2[i]+2*ky3[i]+ky4[i]);  y[i+1]:=y[i]+dy[i]; t:=t+h;  chart1.SeriesList[j].AddXY(t,y[i],'',clblue);  chart2.SeriesList[j].AddXY(t,x[i],'',clblue);    end;    a:=a+da;   j:=j+1;   until a > a1 end; procedure TForm5.SpeedButton1Click(Sender: TObject); begin form4.show;form5.Hide; end;end. unit Unit6; interface 23 var i,p :integer; // переменные, использующиеся в циклах   h,t:real;   // h ­ точность, t ­ шаг   n :integer;  //  число шагов // x0,y0 :real; //  начальные значения x и y   tn, tk : integer; // границы   x,y,kx1,kx2,kx3,kx4,ky1,ky2,ky3,ky4,dx,dy :array[1..10000] of  real;   // x,y ­ значения x и y на соответствующем шаге   // kx1­4,ky1­4 ­ переменные принимающие участие в формуле  Рунге­Кутта   //dx, dy ­ отрезки x и y , из которых складывается решение  уравнения по Рунге­Кутту    a:real;  // a ­ parametr begin tn:=strtoint(edit2.Text);tk:=strtoint(edit3.Text); h:=strtofloat(edit4.Text);x[1]:=strtofloat(edit5.Text); y[1]:=strtofloat(edit6.Text);n:=trunc((tk­tn)/h); t:=tn; chart1.SeriesList[0].Clear;chart2.SeriesList[0].Clear; a:=4;   for i:= 1 to (n) do    begin  kx1[i]:=h*(­a*y[i]­6*x[i]);  kx2[i]:=h*(­a*y[i]­6*(x[i]+kx1[i]/2));  kx3[i]:=h*(­a*y[i]­6*(x[i]+kx2[i]/2));  kx4[i]:=h*(­a*y[i]­6*(x[i]+kx3[i]));  dx[i]:=1/6*(kx1[i]+2*kx2[i]+2*kx3[i]+kx4[i]);  x[i+1]:=x[i]+dx[i]; ky1[i]:=h*x[i];  ky2[i]:=h*(x[i]+ky1[i]/2);  ky3[i]:=h*(x[i]+ky2[i]/2);  ky4[i]:=h*(x[i]+ky3[i]);  dy[i]:=1/6*(ky1[i]+2*ky2[i]+2*ky3[i]+ky4[i]);  y[i+1]:=y[i]+dy[i];  t:=t+h;  chart1.SeriesList[0].AddXY(t,y[i],'',clblue); 24 chart2.SeriesList[0].AddXY(t,x[i],'',clblue);   end;end; 2.3  Лабораторная работа №2       Тема: Методы исследования объектов, динамика которых описывается дифференциальными  уравнениями с использованием программы для моделирования   SciLab ( MatLab)   Цель занятия: 1   Получить   практические   навыки   исследования   систем (объектов), динамика которых описывается дифференциальными уравнениями   2.   Научиться   разрабатывать   алгоритм   и   программу   с   среде моделирования  SciLab ( MatLab)   3. Практически усвоить численные методы Эйлера и Рунге­Кутта для  решения дифференциальных уравнений 1и 2 порядков Задачи  занятия: 1. Разработка алгоритма в виде блок ­ схемы 2. Построение графиков   кривых y(x), dy/dx при параметрах a­ const и var,  3. Анализ результатов исследований. Модели объектов исследования                                   aӳ+bý + cy=f Программа исследования 1. a,c,f ­ const, t – var ( t0 – tk, h= 0.1, 0.01) 2. a,f ­ const, t ­var( t0 – tk, h= 0.1, 0.01), b­ var Пример 1 Условие задачи: Составить алгоритм и проект моделирования объекта, динамика которого описывается   дифференциальным уравнением 1 порядка методом Эйлера 25 ý + 4 ysin(t)­5   программа function yd=f(t, y),yd=5*t­4*y*sin(t),endfunction; y0=5;t0=0;t=0:0.01:3; y=ode(y0,t0,t,f); plot(t,y) Пример 2. Условие задачи: Составить алгоритм и проект моделирования  объекта, динамика которого описывается   дифференциальным    уравнением 1 порядка методом Эйлера   ý + а* ysin(t)­5  , при изменении параметра на 4 значения N=4; //disp(’Vvod  N’); //Цикл для ввода элементов в массиве y. 26 a=1;a1=2; for i=1:N function yd=f(t, y), yd=5*t­a*y*sin(t), endfunction; y0=5;t0=0;t=0:0.01:3; y=ode(y0,t0,t,f); plot(t,y) a=a+a1; end disp(y); Пример 3 Условие задачи: Составить алгоритм и проект моделирования  объекта, динамика которого описывается   дифференциальным    уравнением 2 порядка методом Эйлера ӱ + 0.2*ý+ 4*y=2 27 ion dy=syst(t, y)dy=zeros(2,1); dy(1)=y(2); dy(2)=2­0.2*y(2)­4*y(1) endfunction y0=[0;0];t0=0;t=0:0.1:5; y=ode(y0,t0,t,syst); plot(t,y)     //end //disp(y); 28 Пример 4.  Исследовать объект, динамика которого описывается дифференциальным уравнением                  ӱ + 0.2*ý+ 4*y=2 методом визуального моделирования Исходное уравнение преобразуем к виду Ӱ=2­0.2*ý­4*y 29 30 Пример 5 Условие задачи: Составить алгоритм и проект моделирования  объекта, динамика которого описывается   дифференциальным    уравнением 2 порядка методом Эйлера ӱ + 0.2*ý+ 4*y=2 N=4; //disp(’Vvod  N’); //Цикл для ввода элементов в массиве y. //a=1;a1=5; //for i=1:N function dy=syst(t, y) dy=zeros(2,1); dy(1)=y(2); dy(2)=2­0.2*y(2)­4*y(1) endfunction y0=[0;0];t0=0;t=0:0.1:5; y=ode(y0,t0,t,syst); plot(t,y)     //end 31 //disp(y); Пример. Исследования в среде MatLab Simulink Условие задачи: Составить алгоритм и проект:  моделирования объекта, динамика которого описывается         дифференциальным   уравнением 2 порядка                                      y’’+axy’­y=0,4.   32