|
Рефераты Главная Краткое содержаниепроизведений Архитектура Астрономия Банковское делои кредитование Безопасностьжизнедеятельности Биографии Биология Биржевое дело Бухгалтерия и аудит Военное дело География Геодезия Геология Гражданская оборона Животные Здоровье Земельное право Иностранные языкилингвистика Искусство Историческая личность История История отечественногогосударства и права История политичискихучений История техники Компьютерные сети Компьютеры ЭВМ Криминалистика икриминология Культурология Литература Литература языковедение Маркетинг товароведениереклама Математика Материаловедение Медицина Медицина здоровье отдых Менеджмент (теорияуправления и организации) Металлургия Москвоведение Музыка Наука и техника Нотариат Общениеэтика семья брак Педагогика Право Программированиебазы данных Программное обеспечение Промышленностьсельское хозяйство Психология Радиоэлектроникакомпьютеры и перифирийные устройства Реклама Религия Сексология Социология Теория государства и права Технология Физика Физкультура и спорт Философия Финансовое право Химия - рефераты Хозяйственное право Ценный бумаги Экологическое право Экология Экономика Экономикапредпринимательство Юридическая психология |
Статистический анализ и оптимизация САР Привод сопла ракеты носителяМосковский Государственный Авиационный Институт (Технический университет) Кафедра 704 Информационно-управляющие комплексы Êóðñîâàÿ ðàáîòà Òåìà: “Ñòàòèñòè÷åñêèé àíàëèç ÑÀД Âûïîëíèë: ñò. ãð. 07-403 Êîðíèëîâ Ä.Ì. Ðóêîâîäèòåëü: Êóäðÿøîâ Ñ. Â. Ìîñêâà 1998ã. 1. Задание Ïðèâîä ãèáêîãî ñîïëà ðàêåòû-íîñèòåëÿ: Äàííûå: Ïàðàìåòðû âîçäåéñòâèé íà âõîäàõ ñèñòåìû çàäàíû â âèäå êîððåëÿöèîííûõ ôóíêöèé: ãäå: l=1.2 Óñèëèòåëüíîå çâåíî: Kp=13 1/c Ïàðàìåòðû ïåðâîé íåëèíåéíîñòè: S=0.5; K1=2.5 K2=1.9 Ïàðàìåòðû âòîðîé íåëèíåéíîñòè: S=27 °/c Ïàðàìåòðû òðåòüåé íåëèíåéíîñòè: S=3 ° 1. Ïðîìîäåëèðîâàòü ñîñòîÿíèå ñèñòåìû â íåëèíåéíîì âèäå. 2. Ëèíåàðèçîâàòü íåëèíåéíûå ýëåìåíòû è ïðîìîäåëèðîâàòü ñîñòîÿíèå ñèñòåìû â ëèíåàðèçîâàííîì âèäå. 3. Ïîñòðîèòü ýâîëþöèþ ìàòðèöû êîâàðèàöèé Содержание 1. Задание............................................................................................................................................................. 2. Теоретическая часть............................................................................................................................ 2.1 Случайные процессы и их математическое описание............................................ 2.2 Прохождение стационарного процесса через линейную динамическую систему 2.3 Формирующий фильтр...................................................................................................................... 2.4 Априорный статистический анализ...................................................................................... 2.5 Статистическая линеаризация................................................................................................ 3. Реализация.................................................................................................................................................... 3.1 Система дифференциальных уравнений............................................................................ 3.2 Расчет системы в нелинейной форме.................................................................................... 3.3 Расчет линеаризованной системы......................................................................................... 4. Заключение.................................................................................................................................................. 5. Список литературы................................................................................................................................ 2. Теоретическая часть 2.1 Случайные процессы и их математическое описание Ïóñòü t ïðèíàäëåæèò T (äîïóñòèìîìó ìíîæåñòâó). Åñëè t ïðîáåãàåò íåïðåðûâíûå çíà÷åíèÿ íà ìíîæåñòâå T, òî x(t) ïðèíÿòî íàçûâàòü ñëó÷àéíûì ïðîöåññîì. Ïðè êàæäîì ôèêñèðîâàííîì t=t* âîçíèêàåò ñëó÷àéíàÿ âåëè÷èíà x(t*) êîòîðóþ ïðèíÿòî íàçûâàòü çíà÷åíèåì ñëó÷àéíîãî ïðîöåññà. Ñëó÷àéíûé ïðîöåññ õàðàêòåðèçóåòñÿ ñîâîêóïíîñòüþ ïëîòíîñòåé ðàñïðåäåëåíèÿ âåðîÿòíîñòåé ñ âîçðàñòàþùåé ðàçìåðíîñòüþ k=1,2,...,n. Äåéñòâèòåëüíî âåëè÷èíà ðàâíà âåðîÿòíîñòè òîãî, ÷òî Ïîýòîìó ÷åì áîëüøå n, òåì áîëåå ïîëíîé èíôîðìàöèåé î ïîâåäåíèè x(t) â èíòåðåñóþùåì íàñ èíòåðâàëå âðåìåíè ìû ðàñïîëàãàåì. Ïðàêòè÷åñêè îãðàíè÷èâàþòñÿ ðàññìîòðåíèåì òîëüêî îäíîìåðíûõ è äâóìåðíûõ ïëîòíîñòåé ðàñïðåäåëåíèÿ ëèáî èíûõ õàðàêòåðèñòèê ñëó÷àéíûõ ïðîöåññîâ (ãëàâíûì îáðàçîì ìîìåíòîâ ïåðâîãî è âòîðîãî ïîðÿäêîâ), êîòîðûå îïðåäåëÿþòñÿ äàííûìè ïëîòíîñòÿìè. Ïðèìåðîì ñëó÷àéíîãî ïðîöåññà, ïîëíîñòüþ õàðàêòåðèçóåìîãî îäíîìåðíîé è äâóìåðíîé ïëîòíîñòÿìè, ÿâëÿåòñÿ ìàðêîâñêèé ñëó÷àéíûé ïðîöåññ. Çàâèñèìîñòü ìåæäó çíà÷åíèÿì x(ti) ÿâëÿåòñÿ ïðîñòåéøåé, òàê êàê ðàñïðîñòðàíÿåòñÿ ëèøü íà ñîñåäíèå çíà÷åíèÿ x(ti-1) è x(ti). Íàëè÷èå ïîäîáíîé çàâèñèìîñòè ïðèâîäèò ê òîìó, ÷òî âåðîÿòíîñòü íàõîæäåíèÿ x(ti) â èíòåðâàëå [xi, xi+dxi] â ìîìåíò âðåìåíè t=ti ÿâëÿåòñÿ óñëîâíîé è çàâèñèò îò çíà÷åíèÿ ñëó÷àéíîãî ïðîöåññà â ïðåäûäóùèé ìîìåíò âðåìåíè ti-1.Çàâèñèìîñòü x(ti) îò áîëåå ðàííèõ ìîìåíòîâ âðåìåíè t1, t2, ti-2, (ò. å. îò áîëåå ãëóáîêîé ïðåäûñòîðèè ïðîöåññà) îòñóòñòâóåò. Ýòî îçíà÷àåò, ÷òî äëÿ ìàðêîâñêîãî ïðîöåññà óñëîâíàÿ (èëè ïåðåõîäíàÿ) ïëîòíîñòü: Îòñþäà: Òàêèì îáðàçîì, íà÷àëüíàÿ áåçóñëîâíàÿ îäíîìåðíàÿ ïëîòíîñòü è ñîâîêóïíîñòü óñëîâíûõ (ïåðåõîäíûõ) ïëîòíîñòåé ïîëíîñòüþ îïèñûâàþò ìàðêîâñêèé ñëó÷àéíûé ïðîöåññ. Àáñîëþòíî ñëó÷àéíûì ïðîöåññîì ïðèíÿòî íàçûâàòü òàêîé ïðîöåññ, ëþáûå äâà çíà÷åíèÿ êîòîðîãî ñóòü íåçàâèñèìûå ñëó÷àéíûå âåëè÷èíû.  ýòîì ñëó÷àå ïëîòíîñòü âåðîÿòíîñòè èìååò ñëåäóþùèé âèä: Ñëó÷àéíûé ïðîöåññ íàçûâàåòñÿ ñòàöèîíàðíûì, åñëè âñå åãî ïëîòíîñòè âåðîÿòíîñòåé íå çàâèñÿò îò âûáîðà íà÷àëà îòñ÷åòà âðåìåíè, ò. å. èíâàðèàíòíû ê âðåìåííîìó ñäâèãó t: Èç ýòîãî ñëåäóåò, ÷òî îäíîìåðíàÿ ïëîòíîñòü ðàñïðåäåëåíèÿ ñòàöèîíàðíîãî ïðîöåññà âîîáùå íå çàâèñèò îò âðåìåíè. Ãàóññîâñêèé ïðîöåññ -ýòî òàêîé ñëó÷àéíûé ïðîöåññ ñêîëü óãîäíî ìåðíàÿ ïëîòíîñòü âåðîÿòíîñòè êîòîðîãî ãàóññîâñêàÿ. n-ðàçìåðíîñòü âåêòîðà X, Kx-ìàòðèöà êîâàðèàöèè mx-ìàòåìàòè÷åñêîå îæèëàíèå. Ãàóññîâñêèé ñëó÷àéíûé ïðîöåññ ÿâëÿåòñÿ ñòàöèîíàðíûì è ìàðêîâñêèì. Ê íàèáîëåå âàæíûì ìîìåíòíûì õàðàêòåðèñòèêàì ñòàöèîíàðíîãî ñëó÷àéíîãî ïðîöåññà îòíîñÿòñÿ ìàòåìàòè÷åñêîå îæèäàíèå, äèñïåðñèÿ, êîððåëÿöèîííàÿ ôóíêöèÿ. Ìàòåìàòè÷åñêîå îæèäàíèå õàðàêòåðèçóåò ñðåäíåå òå÷åíèå ïðîöåññà x(t) ïî âðåìåíè. Äèñïåðñèÿ ñëó÷àéíîãî ïðîöåññà Êîððåëÿöèîííàÿ ôóíêöèÿ ãäå , ;ýòà ôóíêöèÿ ïðåäñòàâëÿåò ñîáîé ñðåäíåå ïðîèçâåäåíèå öåíòðèðîâàííûõ çíà÷åíèé ñëó÷àéíîãî ïðîöåññà â ìîìåíòû âðåìåíè t è t+t. Êîððåëÿöèîííàÿ ôóíêöèÿ õàðàêòåðèçóåò ñòåïåíü ëèíåéíîé ñâÿçè (êîððåëÿöèè) ìåæäó çíà÷åíèÿìè ïðîöåññà, îòñòîÿùèìè äðóã îò äðóãà íà âðåìÿ t. Ïðè t=0, êîððåëÿöèîííàÿ ôóíêöèÿ ðàâíà äèñïåðñèè. Ïîíÿòèå êîððåëÿöèîííîé ôóíêöèè ìîæåò áûòü èñïîëüçîâàíî è äëÿ õàðàêòåðèñòèêè ñòåïåíè ñâÿçè äâóõ ñëó÷àéíûõ ïðîöåññîâ x(t) è y(t).  ýòîì ñëó÷àå îíà íàçûâàåòñÿ âçàèìíîé êîððåëÿöèîííîé ôóíêöèåé:  òåîðèè àâòîìàòè÷åñêîãî óïðàâëåíèÿ øèðîêî èñïîëüçóþòñÿ îïèñàíèÿ ñëó÷àéíûõ ïðîöåññîâ â ÷àñòîòíîé îáëàñòè èëè, ïî èíîìó, ñïåêòðàëüíîå ïðåäñòàâëåíèå ñëó÷àéíûõ ïðîöåññîâ. Ðàññìîòðèì ïðåîáðàçîâàíèå Ôóðüå îò êîððåëÿöèîííîé ôóíêöèè: Ïîëó÷åííàÿ ôóíêöèÿ åñòü ÷åòíàÿ âåùåñòâåííàÿ ôóíêöèÿ íàçûâàåìàÿ ñïåêòðàëüíîé ïëîòíîñòüþ ñòàöèîíàðíîãî ñëó÷àéíîãî ïðîöåññà. Ñïðàâåäëèâî îáðàòíîå: Äëÿ ñòàöèîíàðíîãî ñëó÷àéíîãî ïðîöåññà n(t) (ñ íóëåâûì ìàòåìàòè÷åñêèì îæèäàíèåì) òèïà áåëîãî øóìà, êîððåëÿöèîííàÿ ôóíêöèÿ èìååò âèä: ãäå d(t)-äåëüòà-ôóíêöèÿ Äèðàêà, à N -èíòåíñèâíîñòü øóìà. Ñïåêòðàëüíàÿ ïëîòíîñòü ýòîãî ïðîöåññà áóäåò: ÷òî ìîæåò áûòü ïðèíÿòî â êà÷åñòâå îïðåäåëåíèÿ áåëîãî øóìà. Âûðàæåíèå îçíà÷àåò, ÷òî ìîùíîñòü ïàðöèàëüíûõ ñîñòàâëÿþùèõ ñëó÷àéíîãî ïðîöåññà n(t) äëÿ ëþáûõ ÷àñòîò îäíà è òà æå. Ïîýòîìó áåëûé øóì ÿâëÿåòñÿ íàèáîëåå èíòåíñèâíûì âèäîì ïîìåõè. 2.2 Прохождение стационарного процесса через линейную динамическую систему Ðàññìîòðèì ëèíåéíóþ äèíàìè÷åñêóþ ñèñòåìó ñ ïîñòîÿííûìè êîýôôèöèåíòàìè. Íà åå âõîä ïîñòóïàåò ñòàöèîíàðíûé ñëó÷àéíûé ïðîöåññ x(t), íà âûõîäå èìååò ìåñòî ïðîöåññ y(t). Òåîðåòè÷åñêè âûõîäíîé ñëó÷àéíûé ïðîöåññ y(t) ÿâëÿåòñÿ ñòàöèîíàðíûì òîëüêî ïîñëå çàòóõàíèÿ ñâîáîäíûõ êîëåáàíèé â ñèñòåìå, òî åñòü ïðè t®¥. Îäíàêî, â èíæåíåðíûõ ïðèëîæåíèÿõ ìû áóäåì ñ÷èòàòü, ÷òî ïåðåõîäíûé ïðîöåññ â ñèñòåìå çàêàí÷èâàåòñÿ çà âðåìÿ, îïðåäåëÿåìîå â ñîîòâåòñòâèå ñ ïðàâèëàìè òåîðèè óïðàâëåíèÿ. Èíûìè ñëîâàìè ìû áóäåì ñ÷èòàòü ïðîöåññ y(t) ñòàöèîíàðíûì ïî èñòå÷åíèè âðåìåíè çàòóõàíèÿ. Ñïåêòðàëüíàÿ ïëîòíîñòü âûõîäíîãî ïðîöåññà èìååò âèä: Äèñïåðñèÿ âûõîäíîãî ïðîöåññà: Êîððåëÿöèîííàÿ ôóíêöèÿ âûõîäíîãî ïðîöåññà: 2.3 Формирующий фильтр Êàê áûëî ïîêàçàíî âûøå áåëûé øóì èìååò ïîñòîÿííóþ ñïåêòðàëüíóþ ïëîòíîñòü âî âñåì äèàïàçîíå ÷àñòîò. Ñïåêòðàëüíàÿ ïëîòíîñòü âñåõ ôèçè÷åñêè ñóùåñòâóþùèõ ñòàöèîíàðíûõ ñëó÷àéíûõ ïðîöåññîâ ïðåäñòàâëÿþò ñîáîé äðîáíî-ðàöèîíàëüíóþ ôóíêöèþ ÷àñòîòû: Ïðè÷åì ñòåïåíü ïîëèíîìà â çíàìåíàòåëå âûøå ñòåïåíè ïîëèíîìà â ÷èñëèòåëå. Òàêèå ôóíêöèè äîïóñêàþò ôàêòîðèçàöèþ: , ãäå êâàäðàò ìîäóëÿ àìïëèòóäíîé õàðàêòåðèñòèêè íåêîåé ôèêòèâíîé ìèíèìàëüíî-ôàçîâîé äèíàìè÷åñêîé ñèñòåìû, êîòîðóþ â äàëüíåéøåì ìû áóäåì íàçûâàòü ôîðìèðóþùåì ôèëüòðîì ñîîòâåòñòâóþùèì ñïåêòðàëüíîé ïëîòíîñòè íåêîåãî ñëó÷àéíîãî ïðîöåññà. 2.4 Априорный статистический анализ Ïîä àïðèîðíûì ñòàòèñòè÷åñêèì àíàëèçîì (èëè àíàëèçîì òî÷íîñòè) ïîíèìàåòñÿ îïðåäåëåíèå ñòàòèñòè÷åñêèõ õàðàêòåðèñòèê (ìàòåìàòè÷åñêèõ îæèäàíèé, äèñïåðñèé, ñïåêòðàëüíûõ ïëîòíîñòåé, ðàñïðåäåëåíèå âåðîÿòíîñòåé è ò. ï.) êîîðäèíàò óïðàâëÿåìîãî äèíàìè÷åñêîãî îáúåêòà ïî èçâåñòíîìó åãî äèôôåðåíöèàëüíîìó óðàâíåíèþ äâèæåíèÿ è ñòàòèñòè÷åñêèì õàðàêòåðèñòèêàì ñëó÷àéíûõ ôàêòîðîâ. Ïóñòü ëèíåàðèçîâàííûå óðàâíåíèÿ âîçìóùåííîãî äâèæåíèÿ óïðàâëÿåìîãî îáúåêòà èìåþò âèä: ãäå x(t)-âåêòîð ñîñòîÿíèÿ (ôàçîâûé âåêòîð), ðàçìåðíîñòè nx1, A(t)-ìàòðèöà êîýôôèöèåíòîâ, ðàçìåðíîñòè nxn. B(t)- ìàòðèöà êîýôôèöèåíòîâ ,áåëûõ øóìîâ, ðàçìåðíîñòè nxm. n-âåêòîð áåëûõ øóìîâ, ðàçìåðíîñòè mx1. Òîãäà äèôôåðåíöèàëüíûå óðàâíåíèÿ äëÿ âåêòîðà ìàòåìàòè÷åñêèõ îæèäàíèé è ìàòðèöû êîâàðèàöèé èìåþò ñëåäóþùèé âèä: Ðàçìåðíîñòü ìàòðèöû êîâàðèàöèè nxn. N-äèàãîíàëüíàÿ ìàòðèöà èíòåíñèâíîñòåé áåëûõ øóìîâ. Äèôôåðåíöèàëüíûå óðàâíåíèÿ (1)-(3) ðåøàþòñÿ îäíèì èç ÷èñëåííûõ ìåòîäîâ èíòåãðèðîâàíèÿ. Òàêèì îáðàçîì, ìû îïðåäåëÿåì âåêòîð ñîñòîÿíèÿ è ñòàòèñòè÷åñêèå õàðàêòåðèñòèêè ñèñòåìû â ëþáîé ìîìåíò âðåìåíè. Ïåðåä íà÷àëîì èíòåãðèðîâàíèÿ, äîëæíû áûòü èçâåñòíû àïðèîðíûå çíà÷åíèÿ âåêòîðà ñîñòîÿíèÿ, âåêòîðà ìàòåìàòè÷åñêèõ îæèäàíèé è ìàòðèöû êîâàðèàöèé â íà÷àëüíûé ìîìåíò âðåìåíè. 2.5 Статистическая линеаризация Ëåãêî âèäåòü, ÷òî äëÿ ðåøåíèÿ óðàâíåíèé èç ïóíêòà 2.4 íåîáõîäèìû ëèíåéíûå ñèñòåìû óðàâíåíèé. Îäíàêî íà ïðàêòèêå ñèñòåìû óïðàâëåíèÿ ìîãóò ñîäåðæàòü (è ÷àùå âñåãî ñîäåðæàò) íåëèíåéíûå ýëåìåíòû, è óðàâíåíèå äëÿ âåêòîðà ñîñòîÿíèé ïðèíèìàåò âèä:  ýòîì ñëó÷àå ïðèìåíÿåòñÿ ìåòîä ñòàòèñòè÷åñêîé ëèíåàðèçàöèè, êîãäà íåëèíåéíûé ýëåìåíò çàìåíÿåòñÿ ëèíåéíûì â íåêîòîðîì ñìûñëå ýêâèâàëåíòíûì. Ïóñòü íåëèíåéíûé ýëåìåíò èìååò ñëåäóþùèé âèä: Ââåäåì ëèíåéíûé ýëåìåíò ñëåäóþùåãî âèäà: , ãäå Íåîáõîäèìî ÷òîáû âåëè÷èíà íà âûõîäå ëèíåéíîãî ýëåìåíòà áûëà ýêâèâàëåíòíà, â íåêîòîðîì ñìûñëå, âåëè÷èíå íà âûõîäå íåëèíåéíîãî ýëåìåíòà. Ñóùåñòâóþò äâà ïîäõîäà: 1. Êðèòåðèé âèäà: M{z}=M{h} D{z}=D{h} Ôîðìóëû äëÿ êîýôôèöèåíòîâ ñòàòèñòè÷åñêîé ëèíåàðèçàöèè: 2. Âòîðîé ñïîñîá çàêëþ÷àåòñÿ â âûïîëíåíèè êðèòåðèÿ âèäà: M{z}=M{h} D{h-z}®min Êîýôôèöèåíò b âû÷èñëÿåòñÿ ïî ôîðìóëå 3. Реализация Äëÿ ðåøåíèÿ ïîñòàâëåííîé çàäà÷è áûëî íàïèñàíî ïðîãðàììíîå îáåñïå÷åíèå ñ ïîìîùüþ ñðåäû Microsoft Visual C++ 4.0 äëÿ ìàòðè÷íûõ îïåðàöèé, ÷èñëåííûõ ìåòîäîâ èíòåãðèðîâàíèÿ. Îñíîâíàÿ çàäà÷à ðåøàåòñÿ â äâóõ ïðîãðàììàõ, äëÿ ðàñ÷åòà íåëèíåéíîé ñèñòåìû è ëèíåàðèçîâàííîé. 3.1 Система дифференциальных уравнений Äëÿ òîãî ÷òîáû ââåñòè â ñèñòåìó ñëó÷àéíûå âîçìóùåíèÿ ñ òðåáóåìûìè êîððåëÿöèîííûìè ôóíêöèÿìè âîñïîëüçóåìñÿ ïîíÿòèåì ôîðìèðóþùåãî ôèëüòðà, äèíàìè÷åñêîãî çâåíà íà âõîä êîòîðîãî ïîñòóïàåò áåëûé øóì, à íà âûõîäå ïðîöåññ ñ òðåáóåìûìè ïàðàìåòðàìè. Èòàê, ñïåêòðàëüíàÿ ïëîòíîñòü òðåáóåìîãî ïðîöåññà èìååò âèä: Ñîãëàñíî ôîðìóëå ïåðåäàòî÷íàÿ ôóíêöèÿ ôîðìèðóþùåãî ôèëüòðà èìååò âèä: Ïîñëå íåñëîæíûõ ïðåîáðàõîâàíèé ïîëó÷àåì äèôôåðåíöèàëüíûå óðàâíåíèÿ äëÿ ôîðìèðóþùèõ ôèëüòðîâ äâóõ âõîäîâ. Òðåòüåé ôàçîâîé êîîðäèíàòîé áóäåò x, çíà÷åíèå íà âõîäå â òðåòüþ íåëèíåéíîñòü. Åãî óðàâíåíèå èìååò âèä: Ýòè óðàâíåíèÿ è ñîñòàâÿò ñèñòåìó äèôôåðåíöèàëüíûõ óðàâíåíèé: 3.2 Расчет системы в нелинейной форме Äëÿ ðàñ÷åòà äàííîé ñèñòåìû â íåëèíåéíîì âèäå áûëà ðàçðàáîòàíà ñîîòâåòñòâóþùàÿ ïðîãðàììà êîòîðàÿ èíòåãðèðîâàëà ñèñòåìó óðàâíåíèÿ ìåòîäîì Ðóíãå-Êóòòà ÷åòâåðòîãî ïîðÿäêà. Äëÿ ðåøåíèå äèôôåðåíöèàëüíîãî óðàâíåíèÿ ìåòîäîì Ðóíãå-Êóòòà íåîáõîäèìî çíàòü íà÷àëüíûå óñëîâèÿ, òî åñòü çíà÷åíèÿ âåêòîðà ñîñòîÿíèé â íóëåâîé ìîìåíò âðåìåíè, â äàííîì ñëó÷àå îíè íóëåâûå. Ïîñëå îñóùåñòâëåíèÿ èíòåãðèðîâàíèÿ ðåçóëüòàòû áûëè çàïèñàíû â ôàéë, è çàòåì ïîñòðîåíû. Ãðàôèêè èìåþò ñëåäóþùèé âèä: Ðèñóíîê 1. 3.3 Расчет линеаризованной системы Äëÿ ðàñ÷åòà ñèñòåìû â ëèíåàðèçîâàííûì âèäå ïî ôîðìóëàì, íåîáõîäèìî ëèíåàðèçîâàòü ñèñòåìó äèôôåðåíöèàëüíûõ óðàâíåíèé, ò. å. çàìåíèòü òðè íåëèíåéíûõ ýëåìåíòà ëèíåéíûìè êàê ýòî ïîêàçàíî â 2.5. Äëÿ òîãî ÷òîáû ïîëó÷èòü êîýôôèöèåíòû ñòàòèñòè÷åñêîé ëèíåàðèçàöèè íàì íåîáõîäèìî çíàòü ïàðàìåòðû ñëó÷àéíîãî ïðîöåññà íà âõîäå â íåëèíåéíîñòü (ìàòåìàòè÷åñêîå îæèäàíèå è äèñïåðñèþ). Ïàðàìåòð âõîäà â òðåòüþ íåëèíåéíîñòü èìååòñÿ â âåêòîðå ñîñòîÿíèé. Ïîýòîìó åãî äèñïåðñèÿ ïðèñóòñòâóåò â ìàòðèöå êîâàðèàöèè. Äëÿ íàõîæäåíèÿ äèñïåðñèé íà âõîäå â ïåðâóþ è âòîðóþ íåëèíåéíîñòè ïðîèçâåäåì ñëåäóþùèå äåéñòâèÿ. Ïðåîáðàçóåì ñèñòåìó òàêèì îáðàçîì ÷òîáû âõîäîì ñèñòåìû îñòàâàëñÿ îäèí èç âõîäîâ â íàøó ñèñòåìó, à âûõîäîì âõîä â íåëèíåéíîñòü (äâà âõîäà ðàññìàòðèâàþòñÿ îòäåëüíî, à çàòåì èõ äèñïåðñèè ñêëàäûâàþòñÿ). Èòàê, ïîñ÷èòàåì äèñïåðñèè îò âîçäåéñòâèé a è h, íà âõîäàõ ïåðâîé è âòîðîé íåëèíåéíîñòè: Íàçîâåì ñèãíàë íà âõîäå â ïåðâóþ íåëèíåéíîñòü õ, âî âòîðóþ y òîãäà óðàâíåíèå ñâÿçûâàþùèå x è a áóäóò èìåòü âèä: Çäåñü êîýôôèöèåíòû ñòàòèñòè÷åñêîé ëèíåàðèçàöèè ñîîòâåòñòâåííî ïåðâîé è âòîðîé íåëèíåéíîñòåé, äàëåå: Ñëåäîâàòåëüíî ïåðåäàòî÷íàÿ ôóíêöèÿ áóäåò èìåòü âèä: Ñïåêòðàëüíàÿ ïëîòíîñòü ñëó÷àéíîãî ïðîöåññà íà âõîäå a èìååò âèä: Òîãäà äèñïåðñèÿ íà âõîäå ïåðâîé íåëèíåéíîñòè áóäåò èìåòü âèä:  òåîðèè èçâåñòíî, ÷òî èíòåãðàëû âèäà: ãäå:  àíàëèòè÷åñêîé ôîðìå èìåþò âèä: , ãäå Ðàññìîòðèì îòäåëüíî çíàìåíàòåëü íàøåãî èíòåãðàëà, ïðèâåäåííûé ê âèäó: Òîãäà: Îòñþäà ñëåäóåò: Îòêóäà: Ðàññóæäàÿ àíàëîãè÷íûì îáðàçîì ïîëó÷èì îñòàëüíûå äèñïåðñèè: Îòñþäà ñëåäóåò, ÷òî äèñïåðñèÿ íà âõîäå â ïåðâóþ íåëèíåéíîñòü èìååò âèä: Íà âõîäå âòîðîé íåëèíåéíîñòè: Äàëåå íåîáõîäèìî ïîëó÷èòü âûðàæåíèÿ äëÿ ñàìèõ êîýôôèöèåíòîâ ñòàòèñòè÷åñêîé ëèíåàðèçàöèè, âîñïîëüçîâàâøèñü âûâåäåííûìè ðàíüøå ñîîòíîøåíèÿìè. Ëåãêî âèäåòü, ÷òî âñå èíòåãðàëû â ýòèõ ôîðìóëàõ áóäóò èìåòü îäèí èç òðåõ, íèæå ïåðå÷èñëåííûõ, âèäîâ: Ýòè èíòåãðàëû, ñ÷èòàþòñÿ â ÷èñëåííîì âèäå è ïîëó÷àþòñÿ ñ ïîìîùüþ ôóíêöèè îøèáîê è ãàóññîâñêîé ïëîòíîñòè âåðîÿòíîñòè. Îíè ðåàëèçîâàíû â ïðîãðàììå â âèäå ôóíêöèé, òîãäà êîýôôèöèåíòû ñòàòèñòè÷åñêîé ëèíåàðèçàöèè äëÿ ïåðâîãî íåëèíåéíîãî ýëåìåíòà áóäóò èìåòü âèä: K0=k2*J1(0,-s,m,D,-1)+l1*J0(0,-s,m,D,-1)+ k1*J1(-s,s,m,D,0)+k2*J1(s,0,m,D,1)+l2*J0(s,0,m,D,1); K1=(k2*J2(0,-s,m,D,-1)+l1*J1(0,-s,m,D,-1)+k1*J1(-s,s,m,D,0)+ k2*J2(s,0,m,D,1)+l2*J1(s,0,m,D,1)-m*K0)/D; Äëÿ âòîðîãî: K0=-s*J0(0,-s,m,D,-1)+J1(-s,s,m,D,0)+s*J0(s,0,m,D,1); K1=(-s*J1(0,-s,m,D,-1)+J2(-s,s,m,D,0)+s*J1(s,0,m,D,1)+ s*m*J0(0,-s,m,D,-1)-m*J1(-s,s,m,D,0)-s*m*J0(s,0,m,D,1))/D; Äëÿ òðåòüåãî: K0=-s*J0(0,-s,m,D,-1)+J1(-s,s,m,D,0)+s*J0(s,0,m,D,1); K1=(-s*J1(0,-s,m,D,-1)+J2(-s,s,m,D,0)+s*J1(s,0,m,D,1)+ s*m*J0(0,-s,m,D,-1)-m*J1(-s,s,m,D,0)-s*m*J0(s,0,m,D,1))/D; Ëèíåàðèçîâàííàÿ ñèñòåìà äîëæíà èìåòü âèä: Âîñïîëüçîâàâøèñü óðàâíåíèÿìè äëÿ íåëèíåéíîé ñèñòåìû ìàòðèöû ìîæíî çàïèñàòü: Ëåãêî âèäåòü, ÷òî ìàòðèöà A íà êàæäîì øàãå èíòåãðèðîâàíèÿ áóäåò èçìåíÿòüñÿ, â çàâèñèìîñòè îò êîýôôèöèåíòîâ ëèíåàðèçàöèè, êîòîðûå â ñâîþ î÷åðåäü çàâèñÿò îò äèñïåðñèé íà âõîäå íåëèíåéíûõ ýëåìåíòîâ, êîòîðûå çàâèñÿò îò äèñïåðñèé íà âõîäå. Ñëåäîâàòåëüíî íà êàæäîì øàãå èíòåãðèðîâàíèÿ íàøåé ñèñòåìû ìû äîëæíû èíòåãðèðîâàòü äèôôåðåíöèàëüíîå óðàâíåíèå äëÿ ìàòðèöû êîâàðèàöèé K: Ðåàëèçóÿ âñå âûøåñêàçàííîå ïîëó÷èì ãðàôèê èçìåíåíèÿ âåêòîðà ñîñòîÿíèé ëèíåàðèçîâàííîé ñèñòåìû: Ðèñóíîê 2. Ýâîëþöèÿ ìàòðèöû êîâàðèàöèè ïî âõîäíûì âîçäåéñòâèÿì èìååò âèä: Ðèñóíîê 3 Äèñïåðñèÿ âûõîäíîé êîîðäèíàòû: Ðèñóíîê 4 4. Заключение Êàê âèäíî èç ðèñóíêîâ 1 è 2 ñòàòèñòè÷åñêàÿ ëèíåàðèçàöèÿ ïðîâåäåíà ïðàâèëüíî. Ïîâåäåíèå ñèñòåìû â íåëèíåéíîé è ëèíåàðèçîâàííîé ôîðìå ïðèìåðíî îäèíàêîâî. Èç ðèñóíêîâ 3 è 4 ñëåäóåò, ÷òî äèñïåðñèè ñõîäÿòñÿ, ò. å. ÷åðåç íåêîòîðûé ïðîìåæóòîê âðåìåíè ïåðåõîäíûå ïðîöåññû â ñèñòåìå çàêàí÷èâàþòñÿ è íà âûõîäå ñèñòåìû èìååòñÿ ñòàöèîíàðíûé ñëó÷àéíûé ïðîöåññ. 5. Список литературы 1. Êîíñïåêò ëåêöèé 2. Ëåáåäåâ è äð. Ñòàòèñòè÷åñêàÿ äèíàìèêà è îïòèìèçàöèÿ óïðàâëåíèÿ ËÀ. 3. Êîðí Ã. è Ò. Ñïðàâî÷íèê ïî ìàòåìàòèêå 4. Ïîëÿêîâ Þ. Â, Êðóãëîâ È. Þ. Ïðîãðàììèðîâàíèå â ñðåäå Òóðáî Ïàñêàëü 5.5 5. Áàáàê Ñ. Â., Âàñèëüåâ Â. È. è äð. Îñíîâû òåîðèè ìíîãîñâÿçíûõ ñèñòåì àâòîìàòè÷åñêîãî óïðàâëåíèÿ ËÀ. //КуЁёютря Ёрсютр яю ётртшётшчхёъющ фшэрьшъх ЛА //Выяюыэшы туфхэт уЁ. 07-403 //КюЁэшыют Д. М. //Ррёчхт ышэхрЁшчютрээющ ёшётхьы //ЛшэхрЁшчртрээря шётхьр #include <math.h> #include <iostream.h> #include <fstream.h> #include <stdlib.h> #include <time.h> #include "matrix.h" #include "runge.h" #include "random.h" //Зрфрэшх const double Da=0.1; const double Dn=0.03; const double K1=2.5; const double K2=1.9; const double l=1.2; const double Kp=13; const double S1=0.5; const double S2=27; const double S3=3; const double dh=0.01; //Кырёё тхътюЁр ёюётюяэшя схч ъютрЁшрцшш template<class T> class CMyRunge:public CRunge4<T>{ public: CMyRunge(TMatrix<T> t,T h, T x):CRunge4<T>(t,h,x){}; TMatrix<T> f(T,TMatrix<T>); }; //Кырёё ьртЁшцы ъютрЁшрцшш template<class T> class CKovariation:public CRunge4<T>{ public: CKovariation(TMatrix<T> t,T h, T x):CRunge4<T>(t,h,x){}; TMatrix<T> f(T,TMatrix<T>); }; //МртЁшцы уЁртэхэшя TMatrix<double> A(3,3); TMatrix<double> B(3,1); //ВхътюЁ юётюяэшя TMatrix<double> X(3,1); //МртЁшцр ъютрЁшрцшш TMatrix<double> K(3,3); //Сыучрщэыщ тхътюЁ TMatrix<double> s(3,1); //Кюэффшцшхэты трт.ышэхрЁшчрцшш double K11,K21,K31; //ВючтЁрщрхт чэръ x double Sign(double x) { if(x<0) return -1; else return 1; } //ИэтхуЁры юшшсъш double erf(double x) { double x1=fabs(x); const double p=0.3275911; const double a1=0.254829592; const double a2=-0.284496736; const double a3=1.421413741; const double a4=-1.453152027; const double a5=1.061405429; double t=1/(1+p*x1); return Sign(x)*(1-(a1*t+a2*pow(t,2)+a3*pow(t,3)+a4*pow(t,4)+a5*pow(t,5))*exp(-x*x)); } //Пыютэюёть тхЁюятэюётш double KGauss(double x,double m,double D) { return 1/sqrt(2*3.14159*D)*exp(-1/2*(x-m)*(x-m)/D); } //ИэтхуЁрыы //ИэтхуЁры ют p(x) double J0(double Low, double Up, double m, double D, int Inf) { switch( Inf) { case -1:return (1+erf((Up-m)/sqrt(2*D)))/2; case 0: return (erf((Up-m)/sqrt(2*D))-erf((Low-m)/sqrt(2*D)) )/2; case +1:return (1-erf((Low-m)/sqrt(2*D)))/2; } return 0; } //ИэтхуЁры ют x*p(x) double J1(double Low, double Up, double m, double D, int Inf) { switch(Inf){ case -1:return (-D*KGauss(Up,m,D)+m*J0(0,Up,m,D,Inf)); case 0: return(-D*(KGauss(Up,m,D)-KGauss(Low,m,D))+m*J0(Low,Up,m,D,Inf) ); case 1: return(D*KGauss(Low,m,D)+m*J0(Low,0,m,D,Inf)); } return 0; } //ИэтхуЁры ют x*x*p(x) double J2(double Low, double Up, double m, double D, int Inf) { switch(Inf){ case -1:return (-D*(Up*KGauss(Up,m,D)-J0(0,Up,m,D,Inf))+m*J1(0,Up,m,D,Inf)); case 0: { double d; d=-D*(Up*KGauss(Up,m,D)-Low*KGauss(Low,m,D)- J0(Low,Up,m,D,Inf))+m*J1(Low,Up,m,D,Inf); return d;} case 1: return(-D*(-Low*KGauss(Low,m,D)-J0(Low,0,m,D,Inf))+m*J1(Low,0,m,D,Inf)); } return 0; } //Стрт. ышэхрЁшчрцшя //Кюэффшцшхэты яхЁтющ эхышэхщэюётш void Linear1(double m, double D, double& K0,double& K1) { const double s=0.5; const double k1=2.5; const double k2=1.9; double l1=s*(k2-k1); double l2=-l1; K0=k2*J1(0,-s,m,D,-1)+l1*J0(0,-s,m,D,-1)+ k1*J1(-s,s,m,D,0)+k2*J1(s,0,m,D,1)+l2*J0(s,0,m,D,1); K1=(k2*J2(0,-s,m,D,-1)+l1*J1(0,-s,m,D,-1)+k1*J1(-s,s,m,D,0)+ k2*J2(s,0,m,D,1)+l2*J1(s,0,m,D,1)-m*K0)/D; } //Кюэффшцшхэты ттюЁющ эхышэхщэюётш void Linear2(double m, double D, double& K0,double& K1) { const double s=27; K0=-s*J0(0,-s,m,D,-1)+J1(-s,s,m,D,0)+s*J0(s,0,m,D,1); K1=(-s*J1(0,-s,m,D,-1)+J2(-s,s,m,D,0)+s*J1(s,0,m,D,1)+ s*m*J0(0,-s,m,D,-1)-m*J1(-s,s,m,D,0)-s*m*J0(s,0,m,D,1))/D; } //Кюэффшцшхэты тЁхтьхщ эхышэхщэюётш void Linear3(double m, double D, double& K0,double& K1) { const double s=3; K0=-s*J0(0,-s,m,D,-1)+J1(-s,s,m,D,0)+s*J0(s,0,m,D,1); K1=(-s*J1(0,-s,m,D,-1)+J2(-s,s,m,D,0)+s*J1(s,0,m,D,1)+ s*m*J0(0,-s,m,D,-1)-m*J1(-s,s,m,D,0)-s*m*J0(s,0,m,D,1))/D; } //ДшёяхЁёшш эр тхюфх эхышэхщэюётхщ //ПхЁтря //ПрЁрьхтЁы: ъюэффшцшхэты ышэхрЁшчрцшш яхЁтющ ш ттюЁющ эхышэхщэюётхщ double Dx(double k1,double k2) { double d1,d2; d1=1/(2*3.14159)*Kp*Kp*K.Element(1,1)*l/(l+Kp*k1*k2); d2=1/(2*3.14159)*Kp*Kp*K.Element(2,2)/(k1*k2-l)/k1/k2; return (d1+d2); } //ВтюЁря double Dy(double k1,double k2) { double d1,d2; d1=1/(2*3.14159)*k1*k1*Kp*Kp*K.Element(1,1)*l/(l+Kp*k1*k2); d2=1/(3.14159)*K1*Kp*K.Element(2,2)/(l-k1*k2*Kp)/k1/k2; return (d1+d2); } //Рхрышчрцшя ъырёёр тхътюЁр ёюётюяэшя template <class T> TMatrix<T> CMyRunge<T>::f(T h,TMatrix<T> m) { double d=1; TMatrix<T>* t=new TMatrix<T>; s.Element(1,1,B.Element(1,1)*gauss(0,d)); s.Element(2,1,B.Element(2,1)*gauss(0,d)); s.Element(3,1,B.Element(3,1)); *t=A*m+s; return *t; } //Рхрышчрцшя ъырёёр ьртЁшцы ъютрЁшрцшш template <class T> TMatrix<T> CKovariation<T>::f(T h,TMatrix<T> m) { TMatrix<T>* t=new TMatrix<T>(3,3); *t=A*m+m*A.Transpose()+B*B.Transpose(); return *t; } //Пюххрыш main() { //Иэшцшрышчрцшя ьртЁшцы A for(int i=1; i<=3;i++) for(int j=1; j<=3;j++) A.Element(i,j,0); A.Element(1,1,-l); A.Element(2,2,-l); A.Element(3,2,1); //Иэшцшрышчрцшя ьртЁшцы B B.Element(1,1,sqrt(Da*l/3.14159)); B.Element(2,1,sqrt(Dn*l/3.14159)); B.Element(3,1,0); //Иэшцшрышчрцшя ьртЁшцы ъютрЁшрцшш for(i=1; i<=3;i++) for(int j=1; j<=3;j++) if(i!=j)K.Element(i,j,0); K.Element(1,1,Da); K.Element(2,2,Dn); K.Element(3,3,(Da*Kp*Kp+Dn)*10); //Иэшцшрышчрцшя тхътюЁр ёюётюяэшя X for(i=1; i<=3;i++) X.Element(i,1,0); //Кырёё Руэух фыя тхътюЁр ёюётюяэшя CMyRunge<double> Vector(X,dh,0); //Кырёё Руэух фыя ъютрЁшрцшш CKovariation<double> Kov(K,dh,0); //ЗртЁртър double k10,k11,k20,k21; Linear1(0,Da*Kp*Kp,k10,k11); Linear2(0,Da*Kp*Kp,k20,k21); double d1=Dx(k11,k21); double d2=Dy(k11,k21); Linear1(0,d1,k10,k11); Linear2(0,d2,k20,k21); K11=k11; K21=k21; Linear3(0,d2+Dn,k10,k11); K31=k11; A.Element(3,1,K11*K21*Kp); A.Element(3,3,-K11*K21*Kp); int first=1; filebuf f; f.open("x.dat",ios::out ); ostream o(&f); filebuf f1; f1.open("k.dat",ios::out ); ostream o1(&f1); double timec=0; double dt=2.0/100; //Оёэютэющ цшъы for(i=1;i<=200;i++) { timec+=dt; //Зряюыэяхь ьртЁшцу A if (!first){ d1=Dx(K11,K21); d2=Dy(K11,K21); Linear1(0,d1,k10,K11); Linear2(0,d2,k20,K21); Linear3(0,K.Element(3,3),k10,K31); A.Element(3,1,K11*K21*Kp); A.Element(3,3,-K11*K21*Kp); first=0; }; X=Vector.Next(); K=Kov.Next(); o<<X.Element(1,1)<<"\t"<<X.Element(2,1)<<"\t"<< K31*X.Element(3,1)<<"\t"<<timec<<endl; o1<<K.Element(1,1)<<"\t"<<K.Element(2,2)<<"\t"<< K31*K31*K.Element(3,3)<<"\t"<<timec<<endl; } cout<<K; cin>>i; return 0; } //НхышэхрЁшчютрээря шётхьр #include <iostream.h> #include <fstream.h> #include <stdlib.h> #include <time.h> #include "matrix.h" #include "runge.h" #include "random.h" //Зрфрэшх const double Da=0.1; const double Dn=0.03; const double K1=2.5; const double K2=1.9; const double l=1.2; const double Kp=13; const double S1=0.5; const double S2=27; const double S3=3; const double dh=0.001; template<class T> class CMyRunge:public CRunge4<T>{ public: CMyRunge(TMatrix<T> t,T h, T x):CRunge4<T>(t,h,x){}; TMatrix<T> f(T,TMatrix<T>); static T F1(T); static T F2(T); static T F3(T); protected: T k1,k2; }; template <class T> TMatrix<T> CMyRunge<T>::f(T h,TMatrix<T> m) { double d=10; k1=sqrt(Da*l/3.14159); k2=sqrt(Dn*l/3.14159); TMatrix<T>* t=new TMatrix<T>(3,1); t->Element(1,1,-l*m.Element(1,1)+k1*gauss(0,d)); t->Element(2,1,-l*m.Element(2,1)+k2*gauss(0,d)); t->Element(3,1,m.Element(2,1)+F2(F1(Kp*(m.Element(1,1)-m.Element(3,1))))); return *t; } template <class T> T CMyRunge<T>::F1(T x) { if(fabs(x)>=S1) return (K2*x+x/fabs(x)*S1*(K1-K2)); else return(K1*x); return 0; } template <class T> T CMyRunge<T>::F2(T x) { if(x<=-S2) return (-S2); if(fabs(x)<S2) return(x); if(x>=S2) return (S2); return 0; } template <class T> T CMyRunge<T>::F3(T x) { if(x<=-S3) return (-S3); if(fabs(x)<S3) return(x); if(x>=S3) return (S3); return 0; } main() { int i; TMatrix<double> m(3,1); m.Element(1,1,0); m.Element(2,1,0); m.Element(3,1,0); CMyRunge<double> r(m,dh,0); filebuf f; f.open("a.dat",ios::out ); ostream o(&f); double t=0; double dt=2.0/100; for(i=1;i<=200;i++){ t+=dt; m=r.Next(); o<<m.Element(1,1)<<"\t"<<m.Element(2,1)<< "\t"<<CMyRunge<double>::F3(m.Element(3,1))<<"\t"<<t<<endl; //cout<<m.Element(1,1)<<"\t"<<m.Element(2,1)<< //"\t"<<CMyRunge<double>::F3(m.Element(3,1))<<endl; } cin>>i; return 0; } //Matrix Template Library, version 1.00 //(c) by Dmitry Kornilov, 1997 //Designed for Microsoft Visual C++ and Borland C++ #ifndef _IOSTREAM_H #include <iostream.h> #endif #ifndef _MATRIX_H #define _MATRIX_H #endif #define INDEX(i,j,col) (i-1)*col+j-1 template <class T> class TMatrix {public: TMatrix(int ,int ); TMatrix(){col=row=0;buffer=NULL;}; TMatrix(TMatrix<T>&); ~TMatrix(); void Element(int,int,T); T Element(int,int); TMatrix<T> Minor(int, int); TMatrix<T> Invert(); TMatrix<T> Transpose(); T Determinant(); int GetRows(){return row;} int GetCols(){return col;} friend TMatrix operator+(TMatrix<T>&,TMatrix<T>&); friend TMatrix operator+(TMatrix<T>&,T); friend TMatrix operator-(TMatrix<T>&,TMatrix<T>&); friend TMatrix operator*(T,TMatrix<T>&); friend TMatrix operator*(TMatrix<T>&,T); friend TMatrix operator*(const TMatrix<T>& ,const TMatrix<T>& ); friend ostream& operator<< (ostream&, const TMatrix<T>&); void operator=(TMatrix&); protected: T *buffer; int row,col; }; class TMatrixError {public: int Error; char* Message; TMatrixError(int n) {Error=n; }; char* GetString(); }; char* TMatrixError::GetString() { if(Error==5) Message="Index Error"; return Message; } template <class T> TMatrix<T>::TMatrix(int rows,int cols) { row=rows; col=cols; buffer=new T[row*col]; } template <class T> TMatrix<T>::~TMatrix() { delete [] buffer; } template <class T> TMatrix<T>::TMatrix(TMatrix<T> &m) { row=m.row; col=m.col; buffer=new T [row*col]; for(int i=1;i<=row;i++) for(int j=1;j<=col;j++) buffer[(i-1)*col+j-1]=m.buffer[(i-1)*m.col+j-1]; } template <class T>void TMatrix<T>::Element(int i,int j, T f) { buffer[(i-1)*col+j-1]=f; } template <class T>T TMatrix<T>::Element(int i,int j) buffer[(i-1)*col+j-1]; template <class T> TMatrix<T> operator+(TMatrix<T>& m1, TMatrix<T>& m2) { TMatrix<T>& m= *new TMatrix<T>(m1.row,m1.col); for(int i=1;i<=m1.row;i++) for(int j=1;j<=m1.col;j++) m.buffer[(i-1)*m1.col+j-1]=m1.buffer[(i-1)*m1.col+j-1]+ m2.buffer[(i-1)*m1.col+j-1]; //m.Element(i,j,m1.Element(i,j)+m2.Element(i,j)); return m; } template <class T> TMatrix<T> operator+(TMatrix<T>& m1, T f) { TMatrix<T>& m= *new TMatrix<T>(m1.row,m1.col); for(int i=1;i<=m1.row;i++) for(int j=1;j<=m1.col;j++) m.buffer[(i-1)*m1.col+j-1]=m1.buffer[(i-1)*m1.col+j-1]+f; //m.Element(i,j,m1.Element(i,j)+m2.Element(i,j)); return m; } template <class T> TMatrix<T> operator-(TMatrix<T>& m1, TMatrix<T>& m2) { TMatrix<T> *m; m=new TMatrix<T>(m1.row,m1.col); for(int i=1;i<=m1.row;i++) for(int j=1;j<=m1.col;j++) m->buffer[(i-1)*m1.col+j-1]=m1.buffer[(i-1)*m1.col+j-1]- m2.buffer[(i-1)*m1.col+j-1]; return *m; } template <class T> TMatrix<T> operator*(T f, TMatrix<T>& m1) { TMatrix<T> *m; m=new TMatrix<T>(m1.row,m1.col); for(int i=1;i<=m1.row;i++) for(int j=1;j<=m1.col;j++) m->buffer[(i-1)*m1.col+j-1]=m1.buffer[(i-1)*m1.col+j-1]*f; return *m; } template <class T> TMatrix<T> operator*(TMatrix<T>& m1, T f) { TMatrix<T>& m =*new TMatrix<T>(m1.row,m1.col); for(int i=1;i<=m1.row;i++) for(int j=1;j<=m1.col;j++) m.buffer[(i-1)*m1.col+j-1]=m1.buffer[(i-1)*m1.col+j-1]*f; return m; } template <class T> void TMatrix<T>::operator=(TMatrix<T>& m) { row=m.row; col=m.col; if(!buffer) buffer=new T [row*col]; for(int i=1;i<=row;i++) for(int j=1;j<=col;j++) buffer[(i-1)*col+j-1]=m.Element(i,j); } template<class T> TMatrix<T> operator*(const TMatrix<T>& m1,const TMatrix<T>& m2) { int i,j,k; T f3,f2,f1; TMatrix<T>* m =new TMatrix<T>(m1.row,m2.col); for(i=1;i<=m1.row;i++) for(j=1;j<=m2.col;j++) { f3=0; for(k=1;k<=m1.col;k++){ f1=m1.buffer[INDEX(i,k,m1.col)]; f2=m2.buffer[INDEX(k,j,m1.col)]; f3=f3+f1*f2;} m->Element(i,j,f3); } return *m; } template <class T> ostream& operator<< (ostream& o, const TMatrix<T>& m) { for(int i=1;i<=m.row;i++){ for(int j=1;j<=m.col;j++) o<<m.buffer[INDEX(i,j,m.col)]<<"\t"; o<<endl; } o<<endl; return o; } template <class T> TMatrix<T> TMatrix<T>::Minor(int i, int j) { TMatrix<T> *m=new TMatrix<T>(row-1,col-1); int count=0; for (int i1=1;i1<=row;i1++) for (int j1=1;j1<=col;j1++) if ((i1!=i) && (j1!=j)) { m->buffer[count]=buffer[INDEX(i1,j1,col)]; count++; } return *m; } template <class T> T TMatrix<T>::Determinant() { if( col==2) { return buffer[INDEX(1,1,2)]*buffer[INDEX(2,2,2)]- buffer[INDEX(1,2,2)]*buffer[INDEX(2,1,2)];} T f=0; TMatrix<T> p; for(int i=1;i<=col;i++){ p=Minor(i,1); f=f+buffer[INDEX(i,1,col)]*p.Determinant()*sign(i+1); } return f; } int sign(int i) { if( (i % 2)==0) return 1; else return -1; } template <class T> TMatrix<T> TMatrix<T>::Invert() { TMatrix<T> *m=new TMatrix<T>(row,col); T f=Determinant(); for (int i=1;i<=row;i++) for (int j=1;j<=col;j++) m->buffer[INDEX(i,j,col)]=sign(i+j)*(Minor(j,i).Determinant())/f; return *m; } template <class T> TMatrix<T> TMatrix<T>::Transpose() { TMatrix<T> *m=new TMatrix<T>(col,row); for (int i=1;i<=row;i++) for (int j=1;j<=col;j++) m->buffer[INDEX(j,i,row)]=buffer[INDEX(i,j,col)]; return *m; } //Random functions //randomize-чрфрэшх эютющ чртЁртъш //random-ёыучрщэюх чшёыю т фшрярчюэх 0-1 //gauss-эюЁьрыьэюх ЁрёяЁхфхыхэшх ё ьрт юцшфрэшхь // m ш Ё ътрфЁртшчхёъшь ютъыюэхэшхь sigma #ifndef _INC_MATH #include <math.h> #endif #include <stdlib.h> #include <TIME.h> void randomize() { srand( (unsigned)time( NULL ) ); } double random() { double* d=new double; *d=(double)rand()/RAND_MAX; return *d; } double gauss(double m,double sigma) { double a,b,r,sq; do{ a=2*random()-1; b=2*random()-1; r=pow(a,2)+pow(b,2); }while (r>1); sq=sqrt(-2*log(r)/r); return(m+sigma*a*sq); } //Runge-Kutt template //version 1.00 //By Dmitry Kornilov #ifndef _MATRIX_H #include "matrix.h" #endif template <class T> class CRunge4 { public: CRunge4(){}; CRunge4(TMatrix<T>,T,T); virtual TMatrix<T> f(T,TMatrix<T>)=0; TMatrix<T> Next(); private: TMatrix<T> prev; T h; T xprev; }; template <class T> CRunge4<T>::CRunge4(TMatrix<T> m,T f,T x) { h=f; prev=m; xprev=x; } template <class T> TMatrix<T> CRunge4<T>::Next() { TMatrix<T>& m= *new TMatrix<T>(prev.GetRows(),prev.GetCols()); TMatrix<T> k1,k2,k3,k4; k1=f(xprev,prev)*h; k2=f(xprev+h/2,prev+k1*0.5)*h; k3=f(xprev+h/2,prev+k2*0.5)*h; k4=f(xprev+h,prev+k3)*h; m=prev+((T)1/6)*(k1+2.0*k2+2.0*k3+k4); prev=m; xprev+=h; return m; } |
|