Домашнее задание по « Математическое моделирование»

МИНИСТЕРСТВО  ОБРАЗОВАНИЯ И НАУКИ УКРАИНЫ 

ОДЕССКИЙ  НАЦИОНАЛЬНЫЙ МОРСКОЙ УНИВЕРСИТЕТ

Кафедра «Техническая кибернетика» 
 
 
 
 
 
 

Домашнее  задание 

по дисциплине « Математическое моделирование» 
 
 
 
 
 

Выполнил

Ст. 4 к.1 гр. СМФ.

Вересов А.А. 
 
 
 
 
 
 
 

Одесса 2011

СОДЕРЖАНИЕ

С. 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

1  ИДЕНТИФИКАЦИЯ ОБЫКНОВЕННОГО ДИФФЕРЕНЦИАЛЬНОГО  УРАВНЕНИЯ ПЕРВОГО ПОРЯДКА

       1.1 Исходные данные для идентификации

 

Таблица 1.1 –  Исходные данные

i	t	X	Y		i	t	X	Y		i	t	X	Y
1	0	1,0024	0,0227		36	0,7	0,8725	0,90122		71	1,4	1,1813	0,87527
2	0,02	1,0542	0,02982		37	0,72	0,8729	0,87737		72	1,42	1,1826	0,91666
3	0,04	1,0978	0,10894		38	0,74	0,8233	0,92163		73	1,44	1,2062	0,93451
4	0,06	1,1476	0,15046		39	0,76	0,8109	0,92537		74	1,46	1,197	0,91259
5	0,08	1,2001	0,17514		40	0,78	0,7805	0,90337		75	1,48	1,2153	0,93364
6	0,1	1,199	0,18938		41	0,8	0,7719	0,88911		76	1,5	1,1929	0,96321
7	0,12	1,2835	0,23808		42	0,82	0,7136	0,87393		77	1,52	1,2215	0,93485
8	0,14	1,2792	0,27433		43	0,84	0,7468	0,8621		78	1,54	1,2131	0,95966
9	0,16	1,302	0,34212		44	0,86	0,7125	0,86535		79	1,56	1,1898	0,97729
10	0,18	1,3727	0,36766		45	0,88	0,7187	0,84179		80	1,58	1,2181	1,00007
11	0,2	1,377	0,43058		46	0,9	0,714	0,83396		81	1,6	1,2377	0,97675
12	0,22	1,3739	0,44738		47	0,92	0,7081	0,85327		82	1,62	1,2104	1,03089
13	0,24	1,4334	0,50731		48	0,94	0,7177	0,85891		83	1,64	1,173	1,02537
14	0,26	1,3945	0,50186		49	0,96	0,7211	0,82673		84	1,66	1,1735	1,04842
15	0,28	1,421	0,52768		50	0,98	0,7323	0,81018		85	1,68	1,1557	1,009
16	0,3	1,4329	0,57083		51	1	0,6959	0,81768		86	1,7	1,166	1,01218
17	0,32	1,4376	0,62989		52	1,02	0,7336	0,8148		87	1,72	1,1352	1,05393
18	0,34	1,4121	0,65453		53	1,04	0,7361	0,8315		88	1,74	1,1529	1,06542
19	0,36	1,3785	0,6793		54	1,06	0,7597	0,82297		89	1,76	1,1064	1,03543
20	0,38	1,3933	0,68757		55	1,08	0,7609	0,8096		90	1,78	1,1051	1,05661
21	0,4	1,3501	0,70852		56	1,1	0,8162	0,84411		91	1,8	1,0612	1,06439
22	0,42	1,3532	0,78427		57	1,12	0,8569	0,82547		92	1,82	1,0371	1,07777
23	0,44	1,324	0,77278		58	1,14	0,8308	0,82118		93	1,84	1,0319	1,06279
24	0,46	1,2889	0,80012		59	1,16	0,9022	0,84288		94	1,86	1,0007	1,07213
25	0,48	1,2747	0,83881		60	1,18	0,9296	0,81148		95	1,88	0,9879	1,04168
26	0,5	1,2278	0,8668		61	1,2	0,9525	0,84961		96	1,9	0,9801	1,02976
27	0,52	1,1683	0,82894		62	1,22	0,9323	0,85789		97	1,92	0,9683	1,03862
28	0,54	1,1351	0,8563		63	1,24	1,0056	0,81658		98	1,94	0,9302	1,02571
29	0,56	1,1353	0,89112		64	1,26	0,9877	0,86935		99	1,96	0,9637	1,02572
30	0,58	1,0602	0,8551		65	1,28	1,0411	0,82		100	1,98	0,9183	1,0197
31	0,6	1,0414	0,89635		66	1,3	1,0387	0,8449
32	0,62	1,03	0,91656		67	1,32	1,0933	0,8352
33	0,64	0,9941	0,91095		68	1,34	1,0972	0,88124
34	0,66	0,9451	0,89356		69	1,36	1,129	0,87167
35	0,68	0,9414	0,92834		70	1,38	1,1161	0,88537

где   t  - значение безразмерное время ;

    X - значения воздействия (в безразмерном виде);

    Y - значения реакции объекта (в безразмерном виде).

       Данные  таблицы отражены на рис. 1.1.

       1.2  Методика идентификации

 

       Проводится  идентификация линейного обыкновенного  дифференциального уравнения (1.1)

                                                               (1.1)

где  τ - время,

   X(τ) - воздействие,

   Y(τ) - реакция оъекта. 

       Для решения задачи идентификации чаще всего выбирается метод наименьших квадратов с аппроксимацией зависимостей X=f(τ) и Y=f(τ) при котором:

1. проводится аппроксимация зависимостей X=f(t) и Y=f(t) на отрезках оси времени гладкими функциями (полиномы невысоких степеней);
2. для моментов времени путем дифференцирования аппроксимирующих функций определяются производные , .
3. значения функций и производных подставляются в идентифицируемое уравнение и определяется сумма квадратов невязок левой и правой частей уравнения d для всех рассматриваемых моментов времени;
4. значения коэффициентов идентифицируемого дифференциального уравнения определяются путем минимизации суммы квадратов невязок левой и правой частей уравнения.

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

       Для проведения идентификации используется метод аппроксимации на смежных отрезках. Аппроксимация зависимостей X=f(t) и Y=f(t) осуществляется полиномами методом наименьших квадратов.

        

        

                     (1.2)

        

где  τ – независимая переменная;

   i – индекс момента времени на оси основной независимой переменной τ;

Δτ  - отрезок времени на котором проводится аппроксимация;

t = τ-τ₁ - локальная (в пределах nz отрезка) координата времени;

j – индекс момента времени на вспомогательной оси независимой переменной t (в пределах локального отрезка времени);

jm – индекс момента времени конца отрезка (на вспомогательной оси t);

nz – индекс отрезка времени;

   а₀ ¸ а_nf  и b₀ ¸ b_mf  - коэффициенты аппроксимирующих полиномов;

   nf, mf - порядки аппроксимирующих полиномов. 

Выражение для суммы квадратов невязок  по всем рассмотренным зонам имеет  вид:

                                                                                (1.3)

где   m – количество рассмотренных точек всей области определения функции (включая все выделенные отрезки),  

          j – индекс точки.

     Необходимым условием минимума функции δ является равенство нулю ее частных производных:

                                                                                        (1.4)

     Подставив выражение (1.3) в (1.4) можно получить систему  линейных алгебраических уравнений (СЛАУ) (1.5).

                                                             (1.5)

     Решив систему линейных алгебраических уравнений (6) получим значения A, k.

     Проведение идентификации отражено в таблице 1.2.

     Оценка  качества идентификации уравнения (1.1) проводится сравнением заданных значений Y и восстановленных значений Ych. Значения Ych получены при численном решении уравнения (1.1) методом трапеций. 
 

     Аппроксимация на отрезках. 

   

   

      Оценка  качества идентификации приведена  на рис.1.4.

      Численное решение для восстановленной  функции Y=f(t) проводим методом трапеций (1.6).

                                                   (1.6) 

     
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

2 ЧИСЛЕННОЕ РЕШЕНИЕ ОБЫКНОВЕННОГО  ДИФФЕРЕНЦИАЛЬНОГО УРАВНЕНИЯ 1 ПОРЯДКА

 

      Для оценки качества численного решения  уравнения (2.1) используем аналитическое  решение.

                                                               (2.1)

где  τ - время,

   X(τ) - воздействие, 

   Y(τ) - реакция оъекта.

       2.1 Аналитическое решение уравнения  (2.1)

 

      Аналитическое решение уравнения  (2.1) при начальных  условиях  t=0  Y=Y₀ и   воздействии (2.2) имеет вид (2.3).

                                (2.2)

                                                     (2.3)

  где     .

Например, для n=4

  .

       2.2 Численные решение уравнения (2.1)

 

      Рассматривается явная разностная схема.

      

          или                        (2.4)  

     где                                            

     Условие абсолютной устойчивости:  или  .

      Рассматривается аналитико-сеточная схема с постоянным воздействием на отрезке интегрирования.

                                   (2.6)

или   

где   ,   

     Условие абсолютной устойчивости:    выполняется всегда.

     Аналитическое и численное решения приведены  в таблице 2.1. 

Таблица 2.1 – Исходные данные

A=	1	a0=	0,5	C0=	0,34		D1яс=	0,800
k=	1	a1=	0,1	C1=	0,16		D2яс=	0,200
Dt=	0,2	a2=	-0,06	C2=	-0,03		D1ас=	0,818731
Y0=	0	a3=	-0,002	C3=	-0,01		D2ас=	0,090635
		a4=	0,002	C4=	0,002		Ymax=	0,651157
							Ymin=	0

 

Таблица 2.2 - Расчет

i	t	X	Ya	Yяс	Yac	Dяс	Dас	dяс	dас
1	0	0,5	0	0	0	0,0E+00	0,00E+00	0,0%	0,00%
2	0,2	0,5176	0,09235	0,1	0,09223	7,65E-03	-1,26E-04	1,2%	-0,02%
3	0,4	0,5303	0,1707	0,18352	0,17049	1,28E-02	-2,15E-04	2,0%	-0,03%
4	0,6	0,5382	0,2367	0,25288	0,23643	1,62E-02	-2,72E-04	2,5%	-0,04%
5	0,8	0,5414	0,29173	0,30995	0,29142	1,82E-02	-3,03E-04	2,8%	-0,05%
6	1	0,54	0,33692	0,35624	0,33661	1,93E-02	-3,11E-04	3,0%	-0,05%
7	1,2	0,5343	0,37326	0,39299	0,37296	1,97E-02	-3,00E-04	3,0%	-0,05%
8	1,4	0,5246	0,4016	0,42125	0,40133	1,97E-02	-2,74E-04	3,0%	-0,04%
9	1,6	0,5113	0,4227	0,44192	0,42247	1,92E-02	-2,35E-04	3,0%	-0,04%
10	1,8	0,4949	0,43727	0,4558	0,43709	1,85E-02	-1,86E-04	2,8%	-0,03%
11	2	0,476	0,44599	0,46363	0,44586	1,76E-02	-1,29E-04	2,7%	-0,02%
12	2,2	0,4552	0,4495	0,4661	0,44943	1,66E-02	-6,61E-05	2,5%	-0,01%
13	2,4	0,4331	0,44847	0,46391	0,44847	1,54E-02	4,66E-08	2,4%	0,00%
14	2,6	0,4106	0,44358	0,45775	0,44365	1,42E-02	6,79E-05	2,2%	0,01%
15	2,8	0,3886	0,43554	0,44833	0,43567	1,28E-02	1,36E-04	2,0%	0,02%
16	3	0,368	0,42507	0,43639	0,42527	1,13E-02	2,02E-04	1,7%	0,03%
17	3,2	0,3498	0,41298	0,42271	0,41324	9,73E-03	2,65E-04	1,5%	0,04%
18	3,4	0,3351	0,40008	0,40812	0,4004	8,04E-03	3,23E-04	1,2%	0,05%
19	3,6	0,325	0,38727	0,39351	0,38765	6,24E-03	3,74E-04	1,0%	0,06%
20	3,8	0,3209	0,3755	0,37981	0,37592	4,31E-03	4,18E-04	0,7%	0,06%
21	4	0,324	0,36577	0,36803	0,36623	2,25E-03	4,53E-04	0,3%	0,07%
22	4,2	0,3358	0,35916	0,35922	0,35964	5,99E-05	4,77E-04	0,0%	0,07%
23	4,4	0,3577	0,3568	0,35453	0,35729	-2,3E-03	4,89E-04	-0,3%	0,08%
24	4,6	0,3912	0,35991	0,35515	0,3604	-4,8E-03	4,87E-04	-0,7%	0,07%
25	4,8	0,4381	0,36977	0,36237	0,37024	-7,4E-03	4,71E-04	-1,1%	0,07%
26	5	0,5	0,38771	0,37751	0,38815	-1,0E-02	4,39E-04	-1,6%	0,07%
27	5,2	0,5787	0,41517	0,40201	0,41556	-1,3E-02	3,89E-04	-2,0%	0,06%
28	5,4	0,6761	0,45364	0,43735	0,45396	-1,6E-02	3,21E-04	-2,5%	0,05%
29	5,6	0,7941	0,50468	0,4851	0,50491	-2,0E-02	2,33E-04	-3,0%	0,04%
30	5,8	0,9347	0,56995	0,54689	0,57007	-2,3E-02	1,23E-04	-3,5%	0,02%
31	6	1,1	0,65116	0,62445	0,65115	-2,7E-02	-8,72E-06	-4,1%	0,00%