forum.boolean.name - Библиотека для моделирования динамических звеньев

Примеры звеньев - смотри тут.
Если вы незнакомы с Теорией Автоматического Управления, то, возможно, вы несколько не понимаете: о чём идёт речь.
Простым языком теория изложена здесь.
Впрочем, чтобы использовать библиотеку на практике - эти знания вам не нужны.

Сразу пример (к сожалению, вам понадобиться Xors3D). Читаем этот топик.
Вот решение с визуализацией (на скорую руку):

Код:

#include <xors3d.h>

#include "DynUnit.h"



int APIENTRY WinMain(HINSTANCE instance, HINSTANCE prevInstance, LPSTR commandLine, int commandShow)

{

    DynUnitSpace::DynUnit XFilt;

    DynUnitSpace::DynUnit YFilt;



    double T=0.4;//постоянная времени

    double KSI=0.5;//затухание (коэффициент демпфирования)



    double NUM[]={1};

    double DENUM[]={T*T,2*KSI*T,1};



    XFilt.SetParam(NUM,DENUM,sizeof(NUM)/sizeof(double),sizeof(DENUM)/sizeof(double));

    YFilt.SetParam(NUM,DENUM,sizeof(NUM)/sizeof(double),sizeof(DENUM)/sizeof(double));



    xAppTitle("sample");

    int SCRX=800;

    int SCRY=600;

        xGraphics3D(SCRX, SCRY, 32, false, true);

    xMoveMouse(SCRX>>1,SCRY>>1);

    double SEC=xMillisecs()*0.001;

    //"затравочка" - звено не просто усилитель, а с памятью

    XFilt(xMouseX(),SEC);

    YFilt(xMouseY(),SEC);

    int xp,yp;

        while(!xKeyDown(1) || xWinMessage("WM_CLOSE"))

        {

        SEC=xMillisecs()*0.001;

        xp=XFilt(xMouseX(),SEC);

        yp=YFilt(xMouseY(),SEC);



        xCls();

        xOval(xp,yp,20,20,true);

            xRenderWorld();

            xFlip();

        }

    return 0;

}

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

Продублирую код заголовочного файла ниже:

Код:

/*

Пролистав книгу о принципах магии и не взглянув на обложку,

сложно не решить, что это книга о разработке программного обеспечения.

Bruce Tognazzini

*/



/*

НАЗВАНИЕ БИБЛИОТЕКИ:

    DynUnit

    версия 2012.10.04.1

АВТОР:

    impersonalis (icq: 11-999-51-51)

НАЗНАЧЕНИЕ:

    простая библиотека для работы с динамическими звеньями

ЛИНКУЕМЫЕ БИБЛИОТЕКИ:

    -

ДОПОЛНИТЕЛЬНО СМОТРИ:

    -

*/



#ifndef DYNUNIT_H

#define DYNUNIT_H



//пространство библиотеки

namespace DynUnitSpace{



    //функции для работы с матрицами



    //умножение матрицы M (nxn) на столбец Col (n). Результат в R.

    inline void Multip_MatCol(double *(*M),double *Col,int n,double *R);



    //умножение вектора Col (n) на скаляр S. Результат в R.

    inline void Multip_VecScal(double *Col,double S,int n,double *R);



    //умножение строки Row (n) на столбец Col (n). Результат возвращается

    inline double Multip_RowCol(double *Row,double *Col,int n);



    //сложение векторов (n) и (n). Результат в R.

    inline void Add_VecVec(double *V1,double *V2,int n,double *R);



    /*Алгоритм решателя:

    ODE1 - Метод Эйлера        http://ru.wikipedia.org/wiki/Метод_Эйлера

    ODE2 - Heun's method        http://en.wikipedia.org/wiki/Heun's_method

    ODE4-  Метод Рунге — Кутты  http://ru.wikipedia.org/wiki/Метод_Рунге_—_Кутты

    */

    enum ODE_TYPE{ODE1,ODE2,ODE4};



    //Класс "динамическое звено"

    class DynUnit{

        public:

            DynUnit();  //к-тр

            ~DynUnit(); //д-тр



            //Установить параметры звена

            //NUM - коэффициенты числителя (от высшей к низшей степени)

            //DEN - коэффициенты знаменателя (от высшей к низшей степени)

            //NUMs, DENs - размеры соответсвующих векторов

            //Н.у. - нулевые

            //(только для "остановленного" звена)

            //В случае ошибки вернёт false. Если ошибка произошла - состояние звена

            //не теряется (остаётся прежним).

            bool SetParam(double NUM[],

                          double DEN[],

                          int NUMs,

                          int DENs);



            //Установить н.у. ("стартовое" состояние вектора X)

            //(только для "остановленного" звена)

            //В случае ошибки вернёт false. Если ошибка произошла - состояние звена

            //не теряется (остаётся прежним).

            bool Set0(double X0[]);



            //Установить тип решателя (см. описание ODE_TYPE)

            //(только для "остановленного" звена)

            //В случае ошибки вернёт false. Если ошибка произошла - состояние звена

            //не теряется (остаётся прежним).

            bool SetODE(DynUnitSpace::ODE_TYPE OT);



            //"Остановить" звено

            void Stop();



            //Функтор для работы со звеном:

            //U - входное воздействие

            //CurTime - время в сек. (появления этого воздействия)

            //          время не обязательно должно идти от ометки 0.

            //Выход - выход звена

            //          если звено не простой усилитель, то первый такт вернёт 0.

            double operator()(double U,double CurTime);



            //Получить последний выход (тот что вернул функтор при последнем вызове)

            double GetY()const;



        private:



            //очистить внутернние массивы и флаги

            void clear();



            //флаг звено "работает"

            bool RUN;



            //флаг звено - простой усилитель

            bool booster;



            //коэффициент усиления (при booster=true)

            double MULT;



            //вектор состояния

            double *X;



            //числитель

            double *NUM;



            //знаменатель

            double *DEN;



            //старшая степень числителя и знаменателя соответственно

            int m,n;



            //выход звена

            double Y;



            //controllable canonical form

            double *(*A);

            double *B;

            double *C;

            double D;

            //вспомогательные вектора (см. методы)

            double *K1;

            double *K2;

            double *K3;

            double *K4;

            //

            double LastTime;    //"последнее" время

            double LastU;      //"последнее" входное воздействие

            //

            ODE_TYPE ODE;      //тип решателя

            //

            //временные буферы для вычислений

            double *TMP1;

            double *TMP2;

            double *TMP3;

            double *TMP4;

    };

}

#endif // DYNUNIT_H

Динамическое звено описывается:
во frontend: коэффициентами числителя и знаменателя передаточной функции звена.
в backend: пространство состояний.
Для моделирования необходимо реализовать решение диф.ур-ий.
Библиотека предусматривает три метода (чем, точнее, тем дольше одна итерация решения), переключаемых функцией SetODE
Метод Эйлера (дефолтный решатель)
Heun's method
Метод Рунге — Кутты.

Библиотека протестирована, а её выход сравнён с аналогичными реализациями на MATLAB. За это огромное спасибо Felix (реализация, тестирование, советы).

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

Пароль от архива

Скрытый текст (вы должны войти под своим логином или зарегистрироваться и иметь 1 сообщение(ий)):

У вас нет прав, чтобы видеть скрытый текст, содержащийся здесь.