Ошибка регулирования на графике

Вместо введения

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

  • системы стабилизации,
  • системы программного управления,
  • следящие системы

Это достаточно широкий класс систем, которые можно найти где угодно. Но какое это отношение имеет к Unity3D и вероятно к играм в частности? В принципе прямое: в любой игре так или иначе использующей симуляцию как элемент геймплея реализуются САУ, к таким играм относятся, например, Kerbal Space Programm, Digital Combat Simulator (бывший Lock On), Strike Suit Zero и т.д. (кто знает еще примеры — пишите в комментариях). В принципе любая игра, моделирующая реальные физические процессы, в том числе и просто кинематику с динамикой движения, может реализовывать те или иные САУ — этот подход проще, естественнее, а у разработчика уже есть есть набор готовых инструментов, предоставленных всякими Вышнеградскими, Ляпуновыми, Калманами, Чебышевами и прочими Коломогоровами, поэтому можно обойтись без изобретения велосипеда, т.к. его уже изобрели, да так, что получилась отдельная наука: Теория автоматического управления. Главное тут не переусердствовать. Одна тут только проблема: рассказывают про ТАУ не везде, не всем, зачастую мало и не очень понятно.

Немножко теории

Классическая система автоматического управления представленная на следующем рисунке:

image

Ключевым элементом любой САУ является регулятор представляющий из себя устройство, которое следит за состоянием объекта управления и обеспечивает требуемый закон управления. Процесс управления включает в себя: вычисление ошибки управления или сигнала рассогласования e(t) как разницы между желаемой уставкой (set point или SP) и текущей величиной процесса (process value или PV), после чего регулятор вырабатывает управляющие сигналы (manipulated value или MV).

Одной из разновидностью регуляторов является пропорционально-интегрально-дифференцирующий (ПИД) регулятор, который формирует управляющий сигнал, являющийся суммой трёх слагаемых: пропорционального, интегрального и дифференциального.

image

Где, $e(t)$ ошибка рассогласования, а также, $ P = K_p cdot e(t)$ — пропорциональная, $ I = K_i cdot int_0^t e(tau)dtau$ — интегральная, $D = K_d cdot frac{de(t)}{dt}$ — дифференциальная составляющие (термы) закона управления, который в итоговом виде описывается следующими формулами

$ e(t) = SP(t) - PV(t), $

$ MV(t) = underbrace{K_p cdot e(t)}_{P} + underbrace{K_i cdot int_0^t e(tau)dtau}_{I} + underbrace{K_d cdot frac{de(t)}{dt}}_{D}, $

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

Интегральная составляющая I накапливает (интегрирует) ошибку регулирования, что позволяет ПИД-регулятору устранять статическую ошибку (установившуюся ошибку, остаточное рассогласование). Или другими словами: интегральное звено всегда вносит некоторое смещение и если система подвержена некоторыми постоянным ошибкам, то оно их компенсирует (за счет своего смещения). А вот если же этих ошибок нет или они пренебрежительно малы, то эффект будет обратным — интегральная составляющая сама будет вносить ошибку смещения. Именно по этой причине её не используют, например, в задачах сверхточного позиционирования. Ключевым недостатком интегрального закона управления является эффект насыщения интегратора (Integrator windup).

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

Таким образом, в зависимости от ситуации могут применятся П-, ПД-, ПИ- и ПИД-регуляторы, но основным законом управления в основном является пропорциональный (хотя в некоторых специфических задачах и могут использоваться исключительно только звенья дифференциаторов и интеграторов).

Казалось бы, вопрос реализации ПИД-регуляторов уже давно избит и здесь на Хабре есть парочка неплохих статей на эту тему в том числе и на Unity3D, также есть неплохая статья PID Without a PhD (перевод) и цикл статей в журнале «Современные технологии автоматизации» в двух частях: первая и вторая. Также к вашим услугам статья на Википедии (наиболее полную читайте в английском варианте). А на форумах коммьюнити Unity3D нет-нет, да и всплывет PID controller как и на gamedev.stackexchange

При вопрос по реализации ПИД-регуляторов несколько глубже чем и кажется. Настолько, что юных самоделкиных, решивших, реализовать такую схему регулирования ждет немало открытий чудных, а тема актуальная. Так что надеюсь сей опус, кому-нибудь да пригодиться, поэтому приступим.

Попытка номер раз

В качестве примера попытаемся реализовать схему регулирования на примере управления поворотом в простенькой космической 2D-аркаде, по шагам, начиная с самого начала (не забыли, что это туториал?).

Почему не 3D? Потому что реализация не измениться, за исключением того, что придется воротить ПИД-регулятор для контроля тангажа, рысканья и крена. Хотя вопрос корректного применения ПИД-регулирования вместе с кватернионами действительно интересный, возможно в будущем его и освящу, но даже в NASA предпочитают углы Эйлера вместо кватернионов, так что обойдемся простенькой моделью на двухмерной плоскости.

Для начала создадим сам объект игровой объект космического корабля, который будет состоять из собственно самого объекта корабля на верхнем уровне иерархии, прикрепим к нему дочерний объект Engine (чисто спецэффектов ради). Вот как это выглядит у меня:

image

А на сам объект космического корабля накидаем в инспекторе всяческих компонент. Забегая вперед, приведу скрин того, как он будет выглядеть в конце:

image
Но это потом, а пока в нем еще нет никаких скриптов, только стандартный джентльменский набор: Sprite Render, RigidBody2D, Polygon Collider, Audio Source (зачем?).

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

Т.к. помимо самого космического корабля есть куча других, менее умных космических объектов, то сначала опишем родительский класс BaseBody, который в себе будет содержать ссылки на на наши компоненты, методы инициализации и уничтожения, а также ряд дополнительных полей и методов, например для реализации небесной механики:

BaseBody.cs

using UnityEngine;
using System.Collections;
using System.Collections.Generic;

namespace Assets.Scripts.SpaceShooter.Bodies
{
    [RequireComponent(typeof(SpriteRenderer))]
    [RequireComponent(typeof(AudioSource))]
    [RequireComponent(typeof(Rigidbody2D))]
    [RequireComponent(typeof(Collider2D))]

    public class BaseBody : MonoBehaviour
    {
        readonly float _deafultTimeDelay = 0.05f;

[HideInInspector]
        public static List<BaseBody> _bodies = new List<BaseBody>();

        #region RigidBody

        [HideInInspector]
        public Rigidbody2D _rb2d;

        [HideInInspector]
        public Collider2D[] _c2d;

        #endregion

        #region References

        [HideInInspector]
        public Transform _myTransform;

        [HideInInspector]
        public GameObject _myObject;

        /// <summary>
        /// Объект, который появляется при уничтожении
        /// </summary>
        public GameObject _explodePrefab;

        #endregion

        #region  Audio

        public AudioSource _audioSource;

        /// <summary>
        /// Звуки, которые проигрываются при получении повреждения
        /// </summary>
        public AudioClip[] _hitSounds;

        /// <summary>
        /// Звуки, которые проигрываются при появлении объекта
        /// </summary>
        public AudioClip[] _awakeSounds;

        /// <summary>
        /// Звуки, которые воспроизводятся перед смертью
        /// </summary>
        public AudioClip[] _deadSounds;

        #endregion

        #region External Force Variables
        /// <summary>
        /// Внешние силы воздйствующие на объект
        /// </summary>
        [HideInInspector]
        public Vector2 _ExternalForces = new Vector2();

        /// <summary>
        /// Текущий вектор скорости
        /// </summary>
        [HideInInspector]
        public Vector2 _V = new Vector2();

        /// <summary>
        /// Текущий вектор силы гравитации
        /// </summary>
        [HideInInspector]
        public Vector2 _G = new Vector2();
        #endregion

        public virtual void Awake()
        {
            Init();
        }

        public virtual void Start()
        {

        }

        public virtual void Init()
        {
            _myTransform = this.transform;
            _myObject = gameObject;

            _rb2d = GetComponent<Rigidbody2D>();
            _c2d = GetComponentsInChildren<Collider2D>();
            _audioSource = GetComponent<AudioSource>();

            PlayRandomSound(_awakeSounds);

            BaseBody bb = GetComponent<BaseBody>();
            _bodies.Add(bb);
        }

        /// <summary>
        /// Уничтожение персонажа
        /// </summary>
        public virtual void Destroy()
        {
            _bodies.Remove(this);
            for (int i = 0; i < _c2d.Length; i++)
            {
                _c2d[i].enabled = false;
            }
            float _t = PlayRandomSound(_deadSounds);
            StartCoroutine(WaitAndDestroy(_t));
        }

        /// <summary>
        /// Ждем некоторое время перед уничтожением
        /// </summary>
        /// <param name="waitTime">Время ожидания</param>
        /// <returns></returns>
        public IEnumerator WaitAndDestroy(float waitTime)
        {
            yield return new WaitForSeconds(waitTime);

            if (_explodePrefab)
            {
                Instantiate(_explodePrefab, transform.position, Quaternion.identity);
            }

            Destroy(gameObject, _deafultTimeDelay);
        }

        /// <summary>
        /// Проигрывание случайного звука
        /// </summary>
        /// <param name="audioClip">Массив звуков</param>
        /// <returns>Длительность проигрываемого звука</returns>
        public float PlayRandomSound(AudioClip[] audioClip)
        {
            float _t = 0;
            if (audioClip.Length > 0)
            {
                int _i = UnityEngine.Random.Range(0, audioClip.Length - 1);
                AudioClip _audioClip = audioClip[_i];
                _t = _audioClip.length;
                _audioSource.PlayOneShot(_audioClip);
            }
            return _t;
        }

        /// <summary>
        /// Получение урона
        /// </summary>
        /// <param name="damage">Уровень урона</param>
        public virtual void Damage(float damage)
        {
            PlayRandomSound(_hitSounds);
        }

    }
}

Вроде описали все что надо, даже больше чем нужно (в рамках этой статьи). Теперь отнаследуем от него класс корабля Ship, который должен уметь двигаться и поворачивать:

SpaceShip.cs

using UnityEngine;
using System.Collections;
using System.Collections.Generic;

namespace Assets.Scripts.SpaceShooter.Bodies
{
    public class Ship : BaseBody
    {
        public Vector2 _movement = new Vector2();
        public Vector2 _target = new Vector2();
        public float _rotation = 0f;

        public void FixedUpdate()
        {
            float torque = ControlRotate(_rotation);
            Vector2 force = ControlForce(_movement);

            _rb2d.AddTorque(torque);
            _rb2d.AddRelativeForce(force);
        }

        public float ControlRotate(Vector2 rotate)
        {
            float result = 0f;

            return result;
        }

        public Vector2 ControlForce(Vector2 movement)
        {
            Vector2 result = new Vector2();

            return result;

        }
    }
}

Пока в нем нет ничего интересно, на текущий момент это просто класс-заглушка.

Также опишем базовый(абстрактный) класс для всех контроллеров ввода BaseInputController:

BaseInputController.cs

using UnityEngine;
using Assets.Scripts.SpaceShooter.Bodies;

namespace Assets.Scripts.SpaceShooter.InputController
{
    public enum eSpriteRotation
    {
        Rigth = 0,
        Up = -90,
        Left = -180,
        Down = -270
    }

    public abstract class BaseInputController : MonoBehaviour
    {
        public GameObject _agentObject;
        public Ship _agentBody; // Ссылка на компонент логики корабля
        public eSpriteRotation _spriteOrientation = eSpriteRotation.Up; //Это связано с нестандартной 
                                                                           // ориентации спрайта "вверх" вместо "вправо"

        public abstract void ControlRotate(float dt);
        public abstract void ControlForce(float dt);

        public virtual void Start()
        {
            _agentObject = gameObject;
            _agentBody = gameObject.GetComponent<Ship>();
        }

        public virtual void FixedUpdate()
        {
            float dt = Time.fixedDeltaTime;
            ControlRotate(dt);
            ControlForce(dt);
        }

        public virtual void Update()
        {
            //TO DO
        }
    }
}

И наконец, класс контроллера игрока PlayerFigtherInput:

PlayerInput.cs

using UnityEngine;
using Assets.Scripts.SpaceShooter.Bodies;

namespace Assets.Scripts.SpaceShooter.InputController
{
    public class PlayerFigtherInput : BaseInputController
    {
        public override void ControlRotate(float dt)
        {
            // Определяем позицию мыши относительно игрока
            Vector3 worldPos = Input.mousePosition;
            worldPos = Camera.main.ScreenToWorldPoint(worldPos);

            // Сохраняем координаты указателя мыши
            float dx = -this.transform.position.x + worldPos.x;
            float dy = -this.transform.position.y + worldPos.y;

            //Передаем направление
            Vector2 target = new Vector2(dx, dy);
            _agentBody._target = target;

            //Вычисляем поворот в соответствии с нажатием клавиш
            float targetAngle = Mathf.Atan2(dy, dx) * Mathf.Rad2Deg;
            _agentBody._targetAngle = targetAngle + (float)_spriteOrientation;
        }

        public override void ControlForce(float dt)
        {
            //Передаем movement
            _agentBody._movement = Input.GetAxis("Vertical") * Vector2.up 
                + Input.GetAxis("Horizontal") * Vector2.right;
        }
    }
}

Вроде бы закончили, теперь наконец можно перейти к тому, ради чего все это затевалось, т.е. ПИД-регуляторам (не забыли надеюсь?). Его реализация кажется простой до безобразия:

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace Assets.Scripts.Regulator
{
    [System.Serializable] // Этот атрибут необходим для того что бы поля регулятора 
                                   // отображались в инспекторе и сериализовывались
    public class SimplePID
    {
        public float Kp, Ki, Kd;

        private float lastError;
        private float P, I, D;

        public SimplePID()
        {
            Kp = 1f;
            Ki = 0;
            Kd = 0.2f;
        }

        public SimplePID(float pFactor, float iFactor, float dFactor)
        {
            this.Kp = pFactor;
            this.Ki = iFactor;
            this.Kd = dFactor;
        }

        public float Update(float error, float dt)
        {
            P = error;
            I += error * dt;
            D = (error - lastError) / dt;
            lastError = error;

            float CO = P * Kp + I * Ki + D * Kd;

            return CO;
        }
    }
}

Значения коэффициентов по умолчанию возьмем с потолка: это будет тривиальный единичный коэффициент пропорционального закона управления Kp = 1, небольшое значение коэффициента для дифференциального закона управления Kd = 0.2, который должен устранить ожидаемые колебания и нулевое значение для Ki, которое выбрано потому, что в нашей программной модели нет никаких статичных ошибок (но вы всегда можете их внести, а потом героически побороться с помощью интегратора).

Теперь вернемся к нашему классу SpaceShip и попробуем заюзать наше творение в качестве регулятора поворота космического корабля в методе ControlRotate:

 public float ControlRotate(Vector2 rotate)
 {
      float MV = 0f;
      float dt = Time.fixedDeltaTime;

      //Вычисляем ошибку
      float angleError = Mathf.DeltaAngle(_myTransform.eulerAngles.z, targetAngle);

      //Получаем корректирующее ускорение
      MV = _angleController.Update(angleError, dt);

      return MV;
 }

ПИД-регулятор будет осуществлять точное угловое позиционировая космического корабля только за счет крутящего момента. Все честно, физика и САУ, почти как в реальной жизни.

И без этих ваших Quaternion.Lerp

 if (!_rb2d.freezeRotation)
     rb2d.freezeRotation = true;

 float deltaAngle = Mathf.DeltaAngle(_myTransform.eulerAngles.z, targetAngle);
 float T = dt *  Mathf.Abs( _rotationSpeed / deltaAngle);

 // Трансформируем угол в вектор
Quaternion rot = Quaternion.Lerp(
                _myTransform.rotation,
                Quaternion.Euler(new Vector3(0, 0, targetAngle)),
                T);

 // Изменяем поворот объекта
 _myTransform.rotation = rot;

Получившейся исходный код Ship.cs под спойлером

using UnityEngine;
using Assets.Scripts.Regulator;

namespace Assets.Scripts.SpaceShooter.Bodies
{
    public class Ship : BaseBody
    {
        public GameObject _flame;

        public Vector2 _movement = new Vector2();
        public Vector2 _target = new Vector2();

        public float _targetAngle = 0f;
        public float _angle = 0f;

        [Header("PID")]
        public SimplePID _angleController = new SimplePID();

        public void FixedUpdate()
        {
            float torque = ControlRotate(_targetAngle);
            Vector2 force = ControlForce(_movement);

            _rb2d.AddTorque(torque);
            _rb2d.AddRelativeForce(force);
        }

        public float ControlRotate(float rotate)
        {
            float MV = 0f;
            float dt = Time.fixedDeltaTime;

            _angle = _myTransform.eulerAngles.z;

            //Вычисляем ошибку
            float angleError = Mathf.DeltaAngle(_angle, rotate);

            //Получаем корректирующее ускорение
            MV = _angleController.Update(angleError, dt);

            return MV;
        }

        public Vector2 ControlForce(Vector2 movement)
        {
            Vector2 MV = new Vector2();

            //Кусок кода спецэффекта работающего двигателя ради
            if (movement != Vector2.zero)
            {
                if (_flame != null)
                {
                    _flame.SetActive(true);
                }
            }
            else
            {
                if (_flame != null)
                {
                    _flame.SetActive(false);
                }
            }

            MV = movement;
            return MV;
        }
    }
}

Все? Расходимся по домам?

WTF! Что происходит? Почему корабль поворачивается как-то странно? И почему он так резко отскакивает от других объектов? Неужели этот глупый ПИД-регулятор не работает?

Без паники! Давайте попробуем разобраться что происходит.

В момент получения нового значения SP, происходит резкий (ступенчатый) скачок рассогласования ошибки, которая, как мы помним, вычисляется вот так: $e(t) = SP(t) - PV(t), $ соответственно происходит резкий скачок производной ошибки $frac{de(t)}{dt}$, которую мы вычисляем в этой строчке кода:

D = (error - lastError) / dt;

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

Думаю что настал момент построить графики переходного процесса: ступенчатое воздействие от S(t) = 0 в SP(t) = 90 градусов для тела массой в 1 кг, длинной плеча силы в 1 метр и шагом сетки дифференцирования 0.02 с — прям как в нашем примере на Unity3D (на самом деле не совсем, при построении этих графиков не учитывалось, что момент инерции зависит от геометрии твердого тела, поэтому переходный процесс будет немножко другой, но все же достаточно похожий для демонстрации). Все величены на грифике приведены в абсолютных значениях:
image
Хм, что здесь происходит? Куда улетел отклик ПИД-регулятора?

Поздравляю, мы только что столкнулись с таким явлением как «удар» (kick). Очевидно, что в момент времени, когда процесс еще PV = 0, а уставка уже SP = 90, то при численном дифференцировании получим значение производной порядка 4500, которое умножится на Kd=0.2 и сложится с пропорциональным теромом, так что на выходе мы получим значение углового ускорения 990, а это уже форменное надругательство над физической моделью Unity3D (угловые скорости будут достигать 18000 град/с… я думаю это предельное значение угловой скорости для RigidBody2D).

  • Может стоит подобрать коэффициенты ручками, так чтобы скачок был не таким сильным?
  • Нет! Самое лучше чего мы таким образом сможем добиться — небольшая амплитуда скачка производной, однако сам скачок как был так и останется, при этом можно докрутиться до полной неэффективности дифференциальной составляющей.

Впрочем можете поэкспериментировать.

Попытка номер два. Сатурация

Логично, что привод (в нашем случае виртуальные маневровые двигатели SpaceShip), не может отрабатывать сколько угодно большие значения которые может выдать наш безумный регулятор. Так что первое что мы сделаем — сатурируем выход регулятора:

public float ControlRotate(Vector2 rotate, float thrust)
{
    float CO = 0f;
    float MV = 0f;
    float dt = Time.fixedDeltaTime;

    //Вычисляем ошибку
    float angleError = Mathf.DeltaAngle(_myTransform.eulerAngles.z, targetAngle);

    //Получаем корректирующее ускорение
    CO = _angleController.Update(angleError, dt);

    //Сатурируем
    MV = CO;
    if (MV > thrust) MV = thrust;
    if (MV< -thrust) MV = -thrust;

    return MV;
}

А очередной раз переписанный класс Ship полностью выглядит так

namespace Assets.Scripts.SpaceShooter.Bodies
{
    public class Ship : BaseBody
    {
        public GameObject _flame;

        public Vector2 _movement = new Vector2();
        public Vector2 _target = new Vector2();

        public float _targetAngle = 0f;
        public float _angle = 0f;

        public float _thrust = 1f;

        [Header("PID")]
        public SimplePID _angleController = new SimplePID(0.1f,0f,0.05f);

        public void FixedUpdate()
        {
            _torque = ControlRotate(_targetAngle, _thrust);
            _force = ControlForce(_movement);

            _rb2d.AddTorque(_torque);
            _rb2d.AddRelativeForce(_force);
        }

        public float ControlRotate(float targetAngle, float thrust)
        {
            float CO = 0f;
            float MV = 0f;
            float dt = Time.fixedDeltaTime;

            //Вычисляем ошибку
            float angleError = Mathf.DeltaAngle(_myTransform.eulerAngles.z, targetAngle);

            //Получаем корректирующее ускорение
            CO = _angleController.Update(angleError, dt);

            //Сатурируем
            MV = CO;
            if (MV > thrust) MV = thrust;
            if (MV< -thrust) MV = -thrust;

            return MV;
        }

        public Vector2 ControlForce(Vector2 movement)
        {
            Vector2 MV = new Vector2();

            if (movement != Vector2.zero)
            {
                if (_flame != null)
                {
                    _flame.SetActive(true);
                }
            }
            else
            {
                if (_flame != null)
                {
                    _flame.SetActive(false);
                }
            }

            MV = movement * _thrust;

            return MV;
        }

        public void Update()
        {

        }        
    }
}

Итоговая схема нашего САУ тогда станет уже вот такой
image

При этом уже становится понятно, что выход контроллера CO(t) немного не одно и тоже, что управляемая величина процесса MV(t).

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

Запустив игру, мы обнаружим, что космический корабль стал наконец управляемым:

Если построить графики, то можно увидеть, что реакция контроллера стала уже вот такой:
image
Здесь уже используются нормированные величены, углы поделены на значение SP, а выход контроллера отнормирован относительно максимального значения на котором уже происходит сатурация.

Теперь на графике видно наличие ошибки перерегулирования (overshooting) и затухающие колебания. Уменьшая Kp и увеличивая Kd можно добиться уменьшения колебаний, но зато увеличится время реакции контроллера (скорость поворота корабля). И наоборот, увеличивая Kp и уменьшая Kd — можно добиться увеличения скорости реакции контроллера, но появятся паразитные колебания, которые при определенных (критических) значениях, перестанут быть затухающими.

Ниже приведена известна таблица влияния увеличения параметров ПИД-регулятора (как уменьшить шрифт, а то таблица безе переносов не лезет?):

А общий алгоритм ручной настройки ПИД-регулятора следующий:

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

Каких-то общих значений параметров ПИД-регулятора нет: конкретные значения зависят исключительно от параметров процесса (его передаточной характеристики): ПИД-регулятор отлично работающий с одним объектом управления окажется неработоспособным с другим. Более того, коэффициенты при пропорциональной, интегральной и дифференциальной составляющих еще и взаимозависимы.

В общем не будем о грустном, дальше нас ждет самое интересное…

Попытка номер три. Еще раз производные

Приделав костыль в виде ограничения значений выхода контроллера мы так и не решили самую главную проблему нашего регулятора — дифференциальная составляющая плохо себя чувствует при ступенчатом изменении ошибки на входе регуляторе. На самом деле есть множество других костылей, например, в момент скачкообразного изменения SP «отключать» дифференциальную составляющую или же поставить фильтры нижних частот между SP(t) и операцией $SP(t)-PV(t)$ за счет которого будет происходить плавное нарастание ошибки, а можно совсем развернуться и впендюрить самый настоящий фильтр Калмана для сглаживания входных данных. В общем костылей много, и добавить наблюдателя конечно хотелось бы, но не в этот раз.

Поэтому снова вернемся к производной ошибки рассогласования и внимательно на неё посмотрим:

$ frac{de(t)}{dt} = frac{d(SP(t)-PV(t))}{dt} = frac{dSP(t)}{dt} - frac{dPV(t)}{dt}, $

Ничего не заметили? Если хорошенько присмотреться, то можно обнаружить, что вообще-то SP(t), не меняется во времени (за исключением моментов ступенчатого изменения, когда регулятор получает новую команду), т.е. её производная равна нулю:

$ frac{dSP(t)}{dt} = 0, $

тогда

$ frac{de(t)}{dt} = - frac{dPV(t)}{dt}, $

Иными словами, вместо производной ошибки, которая дифференцируема не везде мы можем использовать производную от процесса, который в мире классической механики как правило непрерывен и дифференцируем везде, а схема нашей САУ уже приобретет следующий вид:
image

$ e(t) = SP(t) - PV(t), $

$ CO(t) = underbrace{K_p cdot e(t)}_{P} + underbrace{K_i cdot int_0^t e(tau)dtau}_{I} - underbrace{K_d cdot frac{dPV(t)}{dt}}_{D}, $

Модифицируем код регулятора:

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace Assets.Scripts.Regulator
{
    [System.Serializable]
    public class SimplePID
    {
        public float Kp, Ki, Kd;
        private float P, I, D;

        private float lastPV = 0f;   

        public SimplePID()
        {
            Kp = 1f;
            Ki = 0f;
            Kd = 0.2f;
        }

        public SimplePID(float pFactor, float iFactor, float dFactor)
        {
            this.Kp = pFactor;
            this.Ki = iFactor;
            this.Kd = dFactor;
        }

        public float Update(float error, float PV, float dt)
        {
            P = error;
            I += error * dt;
            D = -(PV - lastPV) / dt;

            lastPV = PV;

            float CO = Kp * P + Ki * I + Kd * D;

            return CO;
        }
    }
}

И немного изменим метод ControlRotate:

public float ControlRotate(Vector2 rotate, float thrust)
{
     float CO = 0f;
     float MV = 0f;
     float dt = Time.fixedDeltaTime;

     //Вычисляем ошибку
     float angleError = Mathf.DeltaAngle(_myTransform.eulerAngles.z, targetAngle);

     //Получаем корректирующее ускорение
     CO = _angleController.Update(angleError, _myTransform.eulerAngles.z, dt);

     //Сатурируем
     MV = CO;
     if (CO > thrust) MV = thrust;
     if (CO < -thrust) MV = -thrust;

     return MV;
}

И-и-и-и… если запустить игру, то обнаружиться, что на самом деле ничего ничего не изменилось с последней попытки, что и требовалось доказать. Однако, если убрать сатурацию, то график реакции регулятора будет выглядеть вот так:
image
Скачок CO(t) по прежнему присутствует, однако он уже не такой большой как был в самом начале, а самое главное — он стал предсказуемым, т.к. обеспечивается исключительно пропорциональной составляющей, и ограничен максимально возможной ошибкой рассогласования и пропорциональным коэффициентом ПИД-регулятора (а это уже намекает на то, что Kp имеет смысл выбрать все же меньше единицы, например, 1/90f), но не зависит от шага сетки дифференцирования (т.е. dt). В общем, я настоятельно рекомендую использовать именно производную процесса, а не ошибки.

Думаю теперь никого не удивит, но таким же макаром можно заменить $K_p cdot e(t)$ на $-K_p cdot PV(t)$, однако останавливаться на этом мы не будем, можете сами поэкспериментировать и рассказать в комментариях, что из этого получилось (самому интересно)

Попытка номер четыре. Альтернативные реализации ПИД-регулятор

Помимо описанного выше идеального представления ПИД-регулятора, на практике часто применяется стандартная форма, без коэффициентов Ki и Kd, вместо которых используются временные постоянные.

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

$ e(t) = SP(t) - PV(t), $

$ CO(t) =CO_{bias} + K_p cdot Bigl(e(t) + frac{1}{T_i} cdot int_0^t e(tau)dtau - T_d cdot frac{dPV(t)}{dt} Bigl), $

где, $T_d= frac{K_d}{K_p}$ — постоянная дифференцирования, влияющая на прогнозирование состояния системы регулятором,
$T_i = frac{K_p}{K_i}$ — постоянная интегрирования, влияющая на интервал усреднения ошибки интегральным звеном.

Основные принципы настройки ПИД-регулятора в стандартной форме аналогичны идеализированному ПИД-регулятору:

  • увеличение пропорционального коэффициента увеличивает быстродействие и снижает запас устойчивости;
  • с уменьшением интегральной составляющей ошибка регулирования с течением времени уменьшается быстрее;
  • уменьшение постоянной интегрирования уменьшает запас устойчивости;
  • увеличение дифференциальной составляющей увеличивает запас устойчивости и быстродействие

Исходный код стандартной формы, вы можете найти под спойлером

namespace Assets.Scripts.Regulator
{
    [System.Serializable]    
    public class StandartPID
    {
        public float Kp, Ti, Td;
        public float error, CO;
        public float P, I, D;

        private float lastPV = 0f;

        public StandartPID()
        {
            Kp = 0.1f;
            Ti = 10000f;
            Td = 0.5f;
            bias = 0f;
        }

        public StandartPID(float Kp, float Ti, float Td)
        {
            this.Kp = Kp;
            this.Ti = Ti;
            this.Td = Td;
        }

        public float Update(float error, float PV, float dt)
        {
            this.error = error;
            P = error;
            I += (1 / Ti) * error * dt;
            D = -Td * (PV - lastPV) / dt;

            CO = Kp * (P + I + D);
            lastPV = PV;

            return CO;
        }
    }
}

В качестве значений по умолчанию, выбраны Kp = 0.01, Ti = 10000, Td = 0.5 — при таких значениях корабль поворачивается достаточно быстро и обладает некоторым запасом устойчивости.

Помимо такой формы ПИД-регулятора, часто используется т.н. реккурентная форма:

$ CO(t_k)=CO(t_{k-1})+K_pleft[left(1+dfrac{Delta t}{T_i}+dfrac{T_d}{Delta t}right)e(t_k)+left(-1-dfrac{2T_d}{Delta t}right)e(t_{k-1})+dfrac{T_d}{Delta t}e(t_{k-2})right] $

Не будем на ней останавливаться, т.к. она актуальна прежде всего для хардверных программистов, работающих с FPGA и микроконтроллерами, где такая реализация значительно удобнее и эффективнее. В нашем же случае — давайте что-нибудь сваям на Unity3D — это просто еще одна реализация ПИД-контроллера, которая ни чем не лучше других и даже менее понятная, так что еще раз дружно порадуемся как хорошо программировать в уютненьком C#, а не в жутком и страшном VHDL, например.

Вместо заключения. Куда бы еще присобачить ПИД-регулятор

Теперь попробуем немного усложнить управление корабля используя двухконтурное управление: один ПИД-регулятор, уже знакомый нам _angleController, отвечает по прежнему за угловое позиционирование, а вот второй — новый, _angularVelocityController — контролирует скорость поворота:

public float ControlRotate(float targetAngle, float thrust)
{
    float CO = 0f;
    float MV = 0f;
    float dt = Time.fixedDeltaTime;

    _angle = _myTransform.eulerAngles.z;

    //Контроллер угла поворота
    float angleError = Mathf.DeltaAngle(_angle, targetAngle);
    float torqueCorrectionForAngle = 
    _angleController.Update(angleError, _angle, dt);

    //Контроллер стабилизации скорости
    float angularVelocityError = -_rb2d.angularVelocity;
    float torqueCorrectionForAngularVelocity = 
        _angularVelocityController.Update(angularVelocityError, -angularVelocityError, dt);

    //Суммарный выход контроллера
    CO = torqueCorrectionForAngle + torqueCorrectionForAngularVelocity;

    //Дискретизируем с шагом 100            
    CO = Mathf.Round(100f * CO) / 100f;

    //Сатурируем
    MV = CO;
    if (CO > thrust) MV = thrust;
    if (CO < -thrust) MV = -thrust;

    return MV;
}

Назначение второго регулятора — гашение избыточных угловых скоростей, за счет изменения крутящего момента — это сродни наличию углового трения, которое мы отключили еще при создании игрового объекта. Такая схема управления [возможно] позволит получить более стабильное поведение корабля, и даже обойтись только пропорциональными коэффициентами управления — второй регулятор будет гасить все колебания, выполняя функцию, аналогичную дифференциальной составляющей первого регулятора.

Помимо этого, добавим новый класс ввода игрока — PlayerInputCorvette, в котором повороты буду осуществляться уже за счет нажатия клавиш «вправо-влево», а целеуказание с помощью мыши мы оставим для чего-нибудь более полезного, например, для управления турелью. Заодно у нас теперь появился такой параметр как _turnRate — отвечающий за скорость/отзывчивость поворота (не понятно только куда его поместить лучше в InputCOntroller или все же Ship).

public class PlayerCorvetteInput : BaseInputController
{
     public float _turnSpeed = 90f;

     public override void ControlRotate()
     {
         // Находим указатель мыши
         Vector3 worldPos = Input.mousePosition;
         worldPos = Camera.main.ScreenToWorldPoint(worldPos);

         // Сохраняем относительные координаты указателя мыши
         float dx = -this.transform.position.x + worldPos.x;
         float dy = -this.transform.position.y + worldPos.y;

         //Передаем направление указателя мыши
         Vector2 target = new Vector2(dx, dy);
         _agentBody._target = target;

         //Вычисляем поворот в соответствии с нажатием клавиш
         _agentBody._rotation -= Input.GetAxis("Horizontal") * _turnSpeed * Time.deltaTime;
    }

    public override void ControlForce()
    {            
         //Передаем movement
         _agentBody._movement = Input.GetAxis("Vertical") * Vector2.up;
    }
}

Также для наглядности накидаем на коленках скрипт для отображения отладочной информации

namespace Assets.Scripts.SpaceShooter.UI
{
    [RequireComponent(typeof(Ship))]
    [RequireComponent(typeof(BaseInputController))]
    public class Debugger : MonoBehaviour
    {
        Ship _ship;
        BaseInputController _controller;
        List<SimplePID> _pids = new List<SimplePID>();
        List<string> _names = new List<string>();

        Vector2 _orientation = new Vector2();

        // Use this for initialization
        void Start()
        {
            _ship = GetComponent<Ship>();
            _controller = GetComponent<BaseInputController>();

            _pids.Add(_ship._angleController);
            _names.Add("Angle controller");

            _pids.Add(_ship._angularVelocityController);
            _names.Add("Angular velocity controller");

        }

        // Update is called once per frame
        void Update()
        {
            DrawDebug();
        }

        Vector3 GetDiretion(eSpriteRotation spriteRotation)
        {
            switch (_controller._spriteOrientation)
            {
                case eSpriteRotation.Rigth:
                    return transform.right;
                case eSpriteRotation.Up:
                    return transform.up;
                case eSpriteRotation.Left:
                    return -transform.right;
                case eSpriteRotation.Down:
                    return -transform.up;
            }
            return Vector3.zero;
        }

        void DrawDebug()
        {
            //Направление поворота
            Vector3 vectorToTarget = transform.position 
                + 5f * new Vector3(-Mathf.Sin(_ship._targetAngle * Mathf.Deg2Rad), 
                    Mathf.Cos(_ship._targetAngle * Mathf.Deg2Rad), 0f);

            // Текущее направление
            Vector3 heading = transform.position + 4f * GetDiretion(_controller._spriteOrientation);

            //Угловое ускорение
            Vector3 torque = heading - transform.right * _ship._Torque;

            Debug.DrawLine(transform.position, vectorToTarget, Color.white);
            Debug.DrawLine(transform.position, heading, Color.green);
            Debug.DrawLine(heading, torque, Color.red);
        }

        void OnGUI()
        {
            float x0 = 10;
            float y0 = 100;

            float dx = 200;
            float dy = 40;

            float SliderKpMax = 1;
            float SliderKpMin = 0;
            float SliderKiMax = .5f;
            float SliderKiMin = -.5f;
            float SliderKdMax = .5f;
            float SliderKdMin = 0;

            int i = 0;
            foreach (SimplePID pid in _pids)
            {
                y0 += 2 * dy;

                GUI.Box(new Rect(25 + x0, 5 + y0, dx, dy), "");

                pid.Kp = GUI.HorizontalSlider(new Rect(25 + x0, 5 + y0, 200, 10), 
                    pid.Kp, 
                    SliderKpMin, 
                    SliderKpMax);
                pid.Ki = GUI.HorizontalSlider(new Rect(25 + x0, 20 + y0, 200, 10), 
                    pid.Ki, 
                    SliderKiMin, 
                    SliderKiMax);
                pid.Kd = GUI.HorizontalSlider(new Rect(25 + x0, 35 + y0, 200, 10), 
                    pid.Kd, 
                    SliderKdMin, 
                    SliderKdMax);

                GUIStyle style1 = new GUIStyle();
                style1.alignment = TextAnchor.MiddleRight;
                style1.fontStyle = FontStyle.Bold;
                style1.normal.textColor = Color.yellow;
                style1.fontSize = 9;

                GUI.Label(new Rect(0 + x0, 5 + y0, 20, 10), "Kp", style1);
                GUI.Label(new Rect(0 + x0, 20 + y0, 20, 10), "Ki", style1);
                GUI.Label(new Rect(0 + x0, 35 + y0, 20, 10), "Kd", style1);

                GUIStyle style2 = new GUIStyle();
                style2.alignment = TextAnchor.MiddleLeft;
                style2.fontStyle = FontStyle.Bold;
                style2.normal.textColor = Color.yellow;
                style2.fontSize = 9;

                GUI.TextField(new Rect(235 + x0, 5 + y0, 60, 10), pid.Kp.ToString(), style2);
                GUI.TextField(new Rect(235 + x0, 20 + y0, 60, 10), pid.Ki.ToString(), style2);
                GUI.TextField(new Rect(235 + x0, 35 + y0, 60, 10), pid.Kd.ToString(), style2);

                GUI.Label(new Rect(0 + x0, -8 + y0, 200, 10), _names[i++], style2);
            }
        }
    }
}

Класс Ship также претерпел необратимые мутации и теперь должен выглядеть вот так:

namespace Assets.Scripts.SpaceShooter.Bodies
{
    public class Ship : BaseBody
    {
        public GameObject _flame;

        public Vector2 _movement = new Vector2();
        public Vector2 _target = new Vector2();

        public float _targetAngle = 0f;
        public float _angle = 0f;

        public float _thrust = 1f;

        [Header("PID")]
        public SimplePID _angleController = new SimplePID(0.1f,0f,0.05f);
        public SimplePID _angularVelocityController = new SimplePID(0f,0f,0f);

        private float _torque = 0f;
        public float _Torque
        {
            get
            {
                return _torque;
            }
        }

        private Vector2 _force = new Vector2();
        public Vector2 _Force
        {
            get
            {
                return _force;
            }
        }

        public void FixedUpdate()
        {
            _torque = ControlRotate(_targetAngle, _thrust);
            _force = ControlForce(_movement, _thrust);

            _rb2d.AddTorque(_torque);
            _rb2d.AddRelativeForce(_force);
        }

        public float ControlRotate(float targetAngle, float thrust)
        {
            float CO = 0f;
            float MV = 0f;
            float dt = Time.fixedDeltaTime;

            _angle = _myTransform.eulerAngles.z;

            //Контроллер угла поворота
            float angleError = Mathf.DeltaAngle(_angle, targetAngle);
            float torqueCorrectionForAngle = 
                _angleController.Update(angleError, _angle, dt);

            //Контроллер стабилизации скорости
            float angularVelocityError = -_rb2d.angularVelocity;
            float torqueCorrectionForAngularVelocity = 
                _angularVelocityController.Update(angularVelocityError, -angularVelocityError, dt);

            //Суммарный выход контроллера
            CO = torqueCorrectionForAngle + torqueCorrectionForAngularVelocity;

            //Дискретизируем с шагом 100            
            CO = Mathf.Round(100f * CO) / 100f;

            //Сатурируем
            MV = CO;
            if (CO > thrust) MV = thrust;
            if (CO < -thrust) MV = -thrust;

            return MV;
        }

        public Vector2 ControlForce(Vector2 movement, float thrust)
        {
            Vector2 MV = new Vector2();

            if (movement != Vector2.zero)
            {
                if (_flame != null)
                {
                    _flame.SetActive(true);
                }
            }
            else
            {
                if (_flame != null)
                {
                    _flame.SetActive(false);
                }
            }

            MV = movement * thrust;

            return MV;
        }

        public void Update()
        {

        }        
    }
}

А вот, собственно заключительное видео того, что должно получиться:

К сожалению получилось охватить не все, что хотелось бы, в частности почти не затронут вопрос настройки ПИД-регулятора и практически не освящена интегральная составляющая — фактически приведен пример только для ПД-регулятора. Собственно изначально планировалось несколько больше примеров (круиз-контроль, вращение турели и компенсация вращательного момента), но статья итак уже разбухла, да и вообще:
image

Немного ссылок

  1. Годная статья на английской вики
  2. PID tutorial
  3. ПИД-регуляторы: вопросы реализации. Часть 1
  4. ПИД-регуляторы: вопросы реализации. Часть 2
  5. PID Without a PhD
  6. PID Without a PhD. Перевод
  7. Derivative Action and PID Control
  8. Control System Lab: PID
  9. ПИД-регулятор своими руками
  10. Корректная реализация разностной схемы ПИД регулятора
  11. Программируем квадрокоптер на Arduino (часть 1)
  12. Виртуальный квадрокоптер на Unity + OpenCV (Часть 1)
  13. Поляков К.Ю. Теория автоматического управления для чайников
  14. PID control system analysis, design, and technology
  15. Aidan O’Dwyer. Handbook of PI and PID Controller Tuning Rules (3rd ed.)
  16. PID process control, a “Cruise Control” example
  17. https://www.mathworks.com/discovery/pid-control.html
  18. http://scilab.ninja/study-modules/scilab-control-engineering-basics/module-4-pid-control/
  19. https://sourceforge.net/p/octave/control/ci/default/tree/inst/optiPID.m

Еще немного ссылок на другие примеры
http://luminaryapps.com/blog/use-a-pid-loop-to-control-unity-game-objects/
http://www.habrador.com/tutorials/pid-controller/3-stabilize-quadcopter/
https://www.gamedev.net/articles/programming/math-and-physics/pid-control-of-physics-bodies-r3885/
https://ksp-kos.github.io/KOS/tutorials/pidloops.html

    1. Виды ошибок регулирования и методы их снижения.

Прямые показатели качества подразделяются
на показатели качества динамического
и установившегося режимов.

Показателями качества динамических
режимов определяются из графика
переходного процесса и основными из
них являются (рис.1.42):

перерегулирование или забросσ,
равный максимуму отклонения значения
переходного процесса относительно
установившегося значения процессаhycm;

— время первой установки t1,
определяемое моментом первого пересечения
графиком переходного процесса
установившегося значенияhycm;

— время переходного процесса tПП,
определяемое момент окончательного
входа графика переходного процесса в
зону допуска, равную±5%от
установившегося значения процессаhycm.

Для всех названных динамических
показателей качества невозможно в общем
случае получить формулы для их расчета.
Это является существенным препятствием
для решения задач анализа и синтеза
САУ.

Показателями качества установившихся
режимов являются ошибки регулирования,
равные абсолютной величине разности
между заданным и фактическим значениями
сигналов САУ и которые в зависимости
от вида входного сигнала САУ подразделяются
на статические (εСТ) и
скоростные ошибки (εСК) и
ошибки (εm)
при отработке гармонического входного
сигнала.

Для
всех названных ошибок регулирования
можно в общем случае получить формулы
их расчета.

Из структурной схемы замкнутой САУ
(рис.1.43) следуют выражения передаточной
функции САУ Wε(p)по ошибке и изображенияε(р)ошибки
регулирования:

Расчет ошибки εmотработки гармонического входного
сигналаx=Xmsinωt
производится по формуле

где
— модуль комплексного числа.

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

(1.54)

Выражение передаточной функции
разомкнутой САУ в общем случае может
быть приведено к виду:

(1.55)

где К– общий коэффициент усиления
разомкнутой САУ:

ν— порядок астатизма САУ, причемνявляется целым неотрицательным
числом.

Для удобства вычислений по формуле
(1.54) подставим в нее выражение WРАЗ(р)из (1.55) и выполним предельный переход:

(1.56)

Статическая ошибка регулирования εСТрассчитывается при постоянном входном
сигналеx(t)=X=const,
а скоростнаяεСК— при
входном сигналеx=Vt,
изменяющемуся во времени с постоянной
скоростьюV=const.
Далее расчеты статической (εСТ)
и скоростной (εСК) ошибок
выполним раздельно.

Расчеты статической ошибки εСт регулирования

Входной сигнал x(t)=X=constи изображением его является.
В соответствии с (1.56) статическую ошибкуεСТследует вычислять по
формуле

(1.57)

1). Пусть в (1.57) значение порядка νастатизма САУ равно нулю:ν=0. Такая
САУ называется статической. Тогда
статическая ошибкаεСТбудет равна

В статической САУ имеется статическая
ошибка εСТ, которую можно
только уменьшить путем увеличения
общего коэффициента усиленияКразомкнутой САУ, но обратить в ноль ее
нельзя.

2). Пусть в (1.57) значение порядка νастатизма САУ равно 1:ν=1. Такая САУ
называется астатической 1-го порядка.
Тогда статическая ошибкаεСТбудет равна

В астатической САУ 1-го порядка статическая
ошибка εСТравна нулю,
т.е САУ является абсолютно точной. Можно
проверить, что при астатизме САУ выше1, статическая ошибка регулирования
всегда будет нулевой.

Соседние файлы в предмете [НЕСОРТИРОВАННОЕ]

  • #
  • #
  • #
  • #
  • #
  • #
  • #
  • #
  • #
  • #
  • #

А

нализ и синтез систем управления

УДК 681.51

ПОВЫШЕНИЕ ТОЧНОСТИ СИСТЕМ С ПИД-РЕГУЛЯТОРАМИ ПРИ ВНЕШНЕМ ВОЗМУЩЕНИИ

А.Г. Александров, Д.А. Хомутов

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

Ключевые слова: ПИД-регулятор, ограниченное возмущение, частотная идентификация.

ВВЕДЕНИЕ

По различным оценкам [1] ПИД-регуляторы [2—4] составляет 90—95 % всех регуляторов, используемых в промышленности. Однако в сложной продукции, такой как самолеты, корабли, электроприводы и др., ПИД -регуляторы используются редко, и соотношение более сложных регуляторов и ПИД-регуляторов обратное. Это связано с тем, что эта продукция работает в условиях интенсивных внешних возмущений (например, порывы ветра для самолетов, волнение моря для кораблей). Для обеспечения малых ошибок регулирования в этих условиях необходимо увеличение коэффициента усиления регулятора, что приводит к необходимости увеличения числа производных ошибки в регуляторе, которые служат для обеспечения устойчивости системы.

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

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

В настоящей работе исследуются системы с ПИ- и ПИД-регуляторами при внешних возмущениях. Получены зависимости ошибок регулирования от уровня внешних возмущений. Если эти ошибки превышают допустимые, то предлагается использовать ПИД -регулятор с производными ошибки регулирования более первой, реализация которого основана на идентификации объекта.

Существо подхода состоит в следующем. Обычно параметры ПИД-регулятора находятся на основе модели объекта, описываемой уравнениями второго порядка и запаздыванием в управлении. Природа запаздывания различна — физическое (транспортное) запаздывание и модельное запаздывание, которое служит для учета неизвестных «малых» постоянных времени. При больших значениях транспортного запаздывания учет дополнительных производных ошибки малоэффективно. Иная ситуация при модельном запаздывании. Здесь точность регулирования может быть значительно увеличена, если известны значения «малых» постоянных времени, и чем больше этих постоянных времени известно, тем большую точность можно обеспечить, определяя по ним коэффициенты регулятора при второй, третьей и далее производных ошибки в ПИД-регуляторе.

Таким образом, для повышения точности ПИД-регулятора необходима идентификация малых постоянных времени. Это затруднено неизвестными внешними возмущениями. Известно несколько методов идентификации в таких условиях: метод инструментальных переменных [5, 6], рандомизированные алгоритмы [7, 8] и метод конечно-частотной идентификации [9]. В этих методах используют идентифицирующее (испытательное) воздействие на объект. Наиболее развит для практического применения метод конечно-частотной идентификации, в соответствии с которым испытательное воздействие представляет собой сумму гармоник с самонастраиваемыми амплитудами и частотами.

Далее исследуются особенности применения этого метода для некоторого класса объектов.

1. СИСТЕМА УПРАВЛЕНИЯ С ПИД-РЕГУЛЯТОРОМ

Рассмотрим систему управления

йтУ(т) + … + У + У = коИ + / < 1 <0, (1)

«(?) = кр

(?) + а0е(0 + |е(т)5т , (2)

е(?) = ,?(?) — У(?), (3)

где у(г) — измеряемый выход объекта (1), и(г) — управление, формируемое регулятором (2), £(<) — задающее воздействие, коэффициенты объекта dj / = 1, …, т, и ко — неизвестные числа, кр, а0 и ах — коэффициенты ПИД-регулятора, /(г) — неизвестное и неизмеряемое, ограниченное внешнее возмущение, е(<) — ошибка регулирования.

Задающее воздействие £(?) является ступенчатой функцией

£(?) = ^о, г > <0, ,?(<) = 0, г = <о. (4)

Число £0 неизвестно, но удовлетворяет неравенству |£0| < £*, в котором £* — заданное число.

Передаточная функция объекта (1) по управлению имеет вид

к„

П ( т> +1) П ( т> +1)

і = 1 і = п + 1

где п — порядок рабочей модели объекта (1), 7] > 0, / = 1, …, п, — постоянные времени рабочей модели объекта — известные числа, 7] > 0, / = п + 1, …, т, — постоянные времени немоделируемой динамики объекта — неизвестные числа.

Постоянные времени Ті, і = 1, упорядочены следующим образом:

т, объекта

І + 1?

і = 1,

т — 1.

(5)

Внешнее возмущение — гармоническая функция с неизвестной амплитудой и частотой из заданного диапазона частот /(г) = акта/г, |а| т /*, г 1 г0 + грег, юа т (0/ < <оь, где юа и — границы диапазона частот внешнего возмущения, /* — заданное число, грег — время регулирования (время затухания переходного процесса, возбужденного задающим воздействием).

Ошибка регулирования в системе е(<) = е^(<) + < 1 г0 + грег, где е^(<) — ошибка регулирования по задающему воздействию, Е/(г) — ошибка регулирования по внешнему возмущению.

Цель управления заключается в том, чтобы ошибка регулирования е(<) удовлетворяла требованию

? 1 ?л + ? ,

0 рег’

(6)

где е* — заданное число.

Пусть цель управления (6) не достигается ни при каких значениях параметров кр, а0 и ах ПИД-регу-лятора (2). Задача состоит в том, чтобы добиться достижения цели (6) с помощью ПИД-регулятора с производными ошибки.

2. СИСТЕМА С ПИД-РЕГУЛЯТОРОМ С ПРОИЗВОДНЫМИ ОШИБКИ

2.1. ПИД-регулятор с производными ошибки

Будем пользоваться следующим дифференциальным уравнением ПИД-регулятора с порядком производной ошибки, равным п — 1:

= к

g, п — 1 п — 1

ип 1 + … + Jgl« + м(?) =

1£П ( ?) + … + А1Е ( ?) + Д0 £( ?) + | е(т) 5т

где д., і = 0, …, п — 1, и Т, і = 1, …, п — 1, — ко-

I ^

эффициенты регулятора.

Передаточная функция ПИД-регулятора с производными ошибки имеет вид

П (+ 1)

п-1

(7)

где ТР., і = 1, …, п, — постоянные времени ПИД-ре-

Р

гулятора с производными ошибки.

*

£

0

п

0

п

п

Далее будем полагать, что коэффициент регу-

лятора кр выбран как

кР 4 Т к ,

п к0

(8)

а его постоянные времени удовлетворяют условиям

Трі = Ті, і = 1, …, п, Tg < Ти.

(9)

здесь опущены, см. работу [9]), р — число гармоник в сигнале, определяемое как

Р = < п>,

(12)

где операция (х) означает округление х вверх до целого.

Выход объекта подается на фильтр Фурье

2.2. Условие достижения цели управления

Утверждение. Пусть параметры регулятора (7) определяются из выражений (8) и (9), тогда найдется такое малое значение 7^, что для достижения цели (6) при условии (5) достаточно, чтобы выполнялось неравенство

4 Тп / * Т1 + Т2 ■

(10)

Отсюда следует, что всегда найдется такой порядок п рабочей модели объекта (1)—(3), при котором ошибка слежения удовлетворяет требованию к точности (6). Таким образом, если вновь выбранный порядок рабочей модели п окажется неудовлетворительным для этого условия, то необходимо увеличивать порядок рабочей модели и строить более сложный регулятор (7) до тех пор, пока это условие не выполнится.

Доказательство этого утверждения приведено в Приложении.

3. ИДЕНТИФИКАЦИЯ ОБЪЕКТА

Для нахождения постоянных времени регулятора, которые выбираются равными постоянным времени объекта (9), необходимо идентифицировать последние. Приведем особенности метода конечно-частотной идентификации [9] для объекта (1), связанные с видом частотных уравнений, используемых в алгоритме конечно-частотной идентификации.

Для идентификации объекта на него подается гармонический испытательный сигнал в виде

р

и(г) = ^ р^тю^, г0 т г < г0 + т, (11)

к = 1

где параметры р1, …, рр — испытательные амплитуды, ю1, … юр — испытательные частоты, т — время идентификации — самонастраиваются в процессе идентификации (алгоритмы самонастроек

а к (т) = — | У(?^ІПЮк?5ґ,

в к (т) = — [ у(?)ео8юк?5?, к = 1, р. (13)

Ркт -1 к

Для оценивания коэффициентов объекта управления (1) решаются следующие частотные уравнения:

5 04) + 1 =

к0

а к + І Р к

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

(14)

где d(^) = dm 5т + dm_ 15т + … + d15 + 1, di, / = 1, …, т, — оценки коэффициентов dj, / = 1, …, т, объекта (1).

4. ПРИМЕР

4.1. Постановка задачи

Рассмотрим объект управления (1) при т = 4:

d4y(4) + d3 у’ + d2 у + d1 у + у = кои + /,

где d1 = 1,248, d2 = 0,2579, d3 = 0,99*10-2, d4 = = 6,4-10—5, ко = 1.

Его передаточная функция имеет вид

Жо(5) =

( Т15 + 1)( Т25 + 1)( Тз5 + 1)( Т^ + 1)

где 71 = 1, Т2 = 0,2, 73 = 0,04, Т4 = 0,008.

Внешнее возмущение /(г) = акто/г с неизвестной амплитудой а, ограниченной числом /* = 5 (а < 5) и частотой Ю/из диапазона [0,1; 100].

На объект в начальный момент времени подается ступенчатое воздействие (4) величиной £о = 30.

Требуется построить ПИД-регулятор, при котором установившаяся ошибка была бы меньше заданной (е* = 0,7).

0

?л + Т

0

0

4.2. Второй порядок рабочей модели объекта (л = 2, Т, и Т2 известны)

Передаточная функция ПИД-регулятора (7) для такого объекта имеет вид

„(4 = * ( 7Р 15 +1) ( V + 1 ), (15)

р (7^ + 1) 5 ’ V ‘

где параметры кр, 7р1, 7р2 и 7^ вычисляем по формулам (8) и (9): кр = 1,25, 7р1 = 1, 7р2 = 0,2 7^ = = 0,17р2 = 0,02.

Максимальное значение ошибки Е/ = 2,23 до-

стигается на частоте 1,1 с 1.

Ошибка по возмущению, оцененная по выражению (10),

4/*72 _ 4 • 5 • 0,2

Т1 + Т2

1 + 0,2

На рис. 1 приведен график выхода объекта, замкнутого указанным ПИ-регулятором при внешнем возмущении /(г) = 58т(1,3г). Ошибка регулирования на этом графике не превышает значения е = 2,22.

4.3. Идентификация третьей постоянной времени

При п = 3 объект (1) примет вид

dзу + d2у + d1 у + у = кои + /(г), г 1 г0, (16)

где dl = 71 + 72 + 73, d2 = 7172 + 7373 + 7273, d3 = 717273, постоянная времени 73 — неизвестна, ее значение необходимо идентифицировать.

Процедура идентификации состоит из следующих операций.

Операция 1. К объекту прикладывается испытательный сигнал (11), где р вычисляем по формуле (12) как р = (1,5) = 2

и(г) = р^тю^ + р28тю2г,

Рис. 1. Выход системы с ПИД-регулятором (15) при задающем воздействии и возмущении

где

Р1 = 2,267, Р2 = 45,16, ®1 = 1/71 = 1,

ю2 = Ь/72 = 25,

где значение Ь = 5 выбирается достаточно большим. Операция 2. Получаем оценки по формулам (13). Операция 3. Оценки коэффициентов уравнения (16) вычисляем по формулам

^ = ко Р 2 ю 1 ( а 1 + Р 1 ) ~ ко Р 1 ю 2 ( а 2 + Iе 2 ) = 1 255 ю1ю2(ю2 — Ю1)(а 1 + Р1)(а2 + Р2)

«2 П2 , -і2 = а1 + Р 1 — к°“ 1 = 0,2618,

2 п 2 П 2,

41(а 1 + р1)

(17)

л д2 л2 л ~2 л2

С = ко Р 2 ю 1 ( а 1 + Р 1 ) — ко Р 1 ю 2 ( а 2 + Р 2 ) = 0 01011

1 2 2 л 2 л 2 л 2 л 2 5

Ю1 ю2(ю2 — Ю1 )(а1 + Р1)(а2 + Р2)

Эти формулы следуют из уравнений (14). Их вывод приведен в Приложении.

Операция 4. Формируем полином (16)

d (5) = d151 + d2 / + с115 + 1 =

= 0,010Ш3 + 0,2618/ + 1,2555 + 1. Находим его корни и представляем его в виде

с (5) = (71 + 1)(725 + 1)(7з5 + 1) =

= (1,0045 + 1)(0,2005 + 1)(0,05035 + 1).

4.4. Третий порядок рабочей модели объекта (л = 3, Г,, Т2 и Т3 известны)

Полученные постоянные времени объекта используем для построения ПИД-регулятора (7)

%(*) -р

= к ( Тр 15 + 1) (Тр25 + 1 ) ( Тр 35 + 1 )

( Т.Д + 1) 5

(18)

где параметры кр, 7тр1, 7р2 и 7р1 вычислены по формулам (8) и (9):

кр = 4,97, 7р1 = 1,004, 7р2 = 0,2,

Тр3 = 0,0503, Т = 0,00503.

Максимальное значение ошибки е = 0,776 достигается на частоте 1,94 с-1.

Ошибка по возмущению, оцененная по выражению (10),

4/*Тз _ 4 • 5 • 0,0503

Т1 + Т2 1 + 0,2

Рис. 2. Выход системы с ПИД-регулятором (18) при задающем воздействии и возмущении

На рис. 2 приведен график выхода объекта, замкнутого указанным ПИД-регулятором. Ошибка регулирования на этом графике не превышает значения е = 0,683.

ЗАКЛЮЧЕНИЕ

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

ПРИЛОЖЕНИЕ

Доказательство утверждения. Выражение для ошибки 8у (?) имеет вид

= ФуЖ (19)

где Фу (я) — передаточная функция системы, связывающая ошибку регулирования с внешним возмущением (передаточная функция ошибки по возмущению). Покажем, что

Ф/(5) ё( я) [ 1 + Жр ( я) Ж0 ( я) ].

Действительно, преобразуя по Лапласу уравнения (1) и (3), а также учитывая передаточную функцию (7), получим следующее выражение для ошибки регулирования:

1

1

1 + Щ( я) ё (^)[ 1 + Щ,( я )Щ,( ^)]

/

Первое слагаемое в правой части этого выражения представляет собой ошибку е^.(0, обусловленную переходным процессом от ступенчатого воздействия (4), второе — ошибку регулирования (?), обусловленную влиянием внешнего возмущения f(?).

Из выражения (19) следует, что цель управления (6) достигается, если выполняется условие

где

н т 6*//*,

Н = Бир |Ф, (»|.

0 < И < ю 7

(20)

(21)

Найдем амплитудно-частотную характеристику |Ф/( У®)1 системы с ПИД-регулятором (7) с производными по ошибке. Нетрудно показать, что

Ф, (я)

= 1 / ё ( я) = _1_

1

7 Щ(я) + 1 ё(я) Щр(я) Щ0(я) + 1

_ (7^ + 1)п—^

А*) ,

где £(я) — характеристический полином замкнутой системы

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

п т

т = я П (7> + 1) П (7> + 1)(7/ + 1)» — 1 +

і = 1 і = п + 1

п

+ П (V + 1)к к = крко.

і = 1

С учетом равенств (9)

Ф/(я) =

(7> + 1)п 1 я

П (7> + 1) і = 1

к + я( 7^5 + 1)п 1 П (+ 1)

і = п + 1

Пренебрегая здесь малыми постоянными времени, получим

ф/ (я) = ——-пт——-.

[к + я] П (7> + 1) і = 1

При условии (8) это выражение может быть представлено в виде

ф/ (я) =

(47> + 1) П (Ті я + 1)

(22)

471

і = 1

Покажем, что малые постоянные времени 7І, і = = п + 1, …, т, и 7. будут мало влиять на запас устойчивости по фазе. Для этого запишем передаточную функцию разомкнутой системы

кр П ( 7 я + 1) = і = 1____________________________

Щ(я) = Щр^Щ/я) = к0

, гр Л ч п — 1 п ІИ

(V + 1) я П (7> + 1) П (7> + 1)

і = 1

к

і = п + 1

я( 7^ + 1 )п 1 П (7> + 1)

і = п + 1

п

т

п

п

6 =

т

Ее частота среза близка к значению к, вычисленному по формуле (8):

®ср = к = І/4 Тп.

Действительно, амплитудно-частотная характеристика разомкнутой системы при таком к близка к 1 на частоте среза шср:

ср

і/,)! =

1/4 Tn

n — 2

1 f T

-4 +1

П

І +1

4Tn^16Tn2 J i = n + W 16Tn2

Фазово-частотная характеристика W(s)

ф(щ) = -2 — (n — l)arctg(Tp) — у arctg(Tiщ), (23)

І = n + І

а условие устойчивости с запасом по фазе имеет вид ф(шср) 1 —п + фз, где фз — запас устойчивости системы (1)—(3) по фазе.

Положим фз = п/4. Подставляя это значение в предыдущее выражение и решая его совместно с выражением (23), получим

т

(п — ^аг^Г^ + X ^Т^ < п .

І = n + І

При условии (5)

X аг^(7Хр) < 8 . (24)

I = П + 1

Действительно, ввиду малости углов запишем это выражение в виде

—У Т. <

4Tn . І i 8

n І = n + І

(25)

Сумма Т, г = п + 1, …, т, при условии (5) удовлетворяет (по свойству суммы геометрической прогрессии) неравенству

х Т < 2Тп + 1 < Тпп + 1

Используя это выражение в формуле (25), получаем, что неравенство (24) выполняется. Следовательно, постоянные времени Т, і = п + 1, …, т, и Т& будут мало влиять на запас устойчивости по фазе.

Оценим уровень амплитудно-частотной характеристики (22). Имеем

|Ф/ (>)! =

4Т„щ

л/16Т2щ2 ПІ Ті2Щ2 + 1

(26)

Оценка сверху (мажоранта) этой функции при п 1 2 имеет вид

|ф*(/щ)1 =

4Т„щ

Iгр2 2 7 [7р2 2 7

^11 щ + 1 Т2щ + 1

О m щ m го. (27)

Действительно, разность мажоранты (27) и амплитудно-частотной характеристики (26) при п > 2 всегда положительна (случай п = 2 рассмотрен далее). Нетрудно видеть, что

!®*(>)! — |ф/ (>)! = 74Т2

4Tn щ

22 , щ +

1 nJ Т2щ2 +1

І=З

■Лгт^ п V т2«>2 *1

‘2-2 п. т22

г = 1

Выражение в квадратных скобках положительно при 0 < ш < го, поэтому все выражение положительно. При ш = 0 и ш ^ го оно обращается в нуль. Таким образом, мажоранта выполняется для случая п > 2.

Покажем, что мажоранта (27) справедлива и для случая п = 2. Имеем

|Ф*(/ш)| — |Ф/(»|л = 2 =

4 І2 щ [ Jie: Г22 щ 2 + 1 — л /І22 щ 2 + 1 і

J Ті2 Щ 2 + 1 ^ 2 8 + 6 2 8 +

Отсюда видно, что выражение в квадратных скобках положительно при 0 < ш < го, а при ш = 0 и ш ^ го оно обращается в нуль. Таким образом, мажоранта (27) выполняется и для случая п = 2.

Найдем максимум мажоранты (27). Ее производная

д| Ф * (/ш ) | = 4ТП ( 1 — Т1 Т2 ш )

,/(Т12ш2 + 1 )^(Т22ш2 + 1)3′

Она обращается в нуль при частоте ш = ш* = 1Д/ Т1Т2, и соответственно

sup |Ф*(»| =

О <а <ю

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

4 Tn

Ті + І2 •

Таким образом, максимум амплитудно-частотной характеристики (26) H = sup |Ф (/щ)| удовлетворяет

О <<а<со

условию

4T

н m———n—, О m щ m го, n і 2. (28)

T1 + T2

Максимумы (28) и (21) связаны выражением H =

І = 1

= Н + £, где £ — малое число, зависящее от малых постоянных времени Т. и Г, г = п + 1, т.

При Т. = 0 и Т’ = 0, г = п + 1, т, число £ = 0. По непрерывности заключаем, что всегда существуют достаточно малые Т и Т, г = п + 1, т, такие, что для лю-

2 Л 2 m I ‘т2

1

m

m

m

m

бого наперед заданного сколь угодно малого £* будет выполняться условие £ < £*. Поэтому для проверки условия (20) будем пользоваться выражением е* > Н/ *. Подставляя сюда условие (20), получаем выражение (10).

Вывод формул (17). Уравнения (14) при п = 3 принимают вид

-М — ^2 + М Шк = ко/( «к + 1Рк ) — 1

Отсюда, выделяя мнимую и действительную части,

получаем следующие уравнения

— ё3 ^ Ш! = к0 -2^, (29)

а 1 + р 1

— ^ ш1 = ко — 1, (30)

а 1 + Р1

— Ъ ш + * ш2 = к0 (31)

«2 + в 2

Оценка с12 входит только в уравнение (30). Для определения оценок с11 и ё3 запишем уравнения (29) и (31) в виде следующей системы:

-ko Р 1

— 2 — 2

3 -Si Si d 3 = — 2 P 2 -1 + P1

3 _-S2 s2 d i -ko P 2

-2 -2 — 2 + р2

Решая их, а также выражая из уравнения (30) оценку с12, получаем формулы (17).

ЛИТЕРАТУРА

1. Денисенко В.В. Разновидности ПИД-регуляторов // Автоматизация в промышленности. — 2GG7. — № б. — С. 45—5G.

2. Ротач В.Я. Расчет промышленных автоматических систем регулирования. — М.: Энергия, 1973.

3. Astrom K.J., and Hagglund T. Advanced PID control. — ISA, 2GG6.

4. Voda A.A., and Landau I.D. A method for the auto-calibration of PID-Controllers // Automatic. — 1995. — Vol. 31, N 1. — P. 41—53.

5. Льюнг Л. Идентификация систем. Теория для пользователя. — М.: Наука, 1991.

6. Wong K.Y., Polak E. Identification of linear discrete-time systems using instrumental variable method // IEEE Trans. Automat. Control. — 19б7. — Vol. AC-12. — P. 7G7—718.

7. Граничин О.Н., Поляк Б.Т. Рандомизированные алгоритмы оценивания и оптимизации при почти производных помехах. — М.: Наука, 2GG3.

8. Бунич А.Л., Бахтадзе Н.Н. Синтез и применение дискретных систем управления с идентификатором. — М.: Наука, 2GG3. — 232 с.

9. Александров А.Г. Конечно-частотная идентификация: самонастройка испытательных частот // Робастное управление и частотная идентификация: Сб. науч. тр. — Электросталь, ЭПИ МИСиС, 2GG4. — С. б7—97.

Статья представлена к публикации членом редколлегии

В.Ю. Рутковским.

Александров Альберт Георгиевич — д-р физ.-мат. наук,

вед. науч. сотрудник, И alex7@ipu.rssi.ru,

Хомутов Дмитрий Алексеевич — вед. инженер,

И serrew@yandex.ru,

Институт проблем управления им. В.А. Трапезникова РАН,

г. Москва, в (495) 334-7б-41.

Содержание специального выпуска сборника «Управление большими системами», посвященного математической теории игр и ее приложениям (http://ubs.mtas.ru)

^ Васильев В.А. Об одной аксиоматизации обобщенного расширения Оуэна. — С. 5—17.

^ Винниченко С.В. Непрерывная игра НИМ. — С. 18—31.

^ Кацев И.В., Яновская Е.Б. Промежуточные между пред к- и пред п-ядрами решения кооперативных игр. — С. 32—54.

^ Мазалов В.В., Сакагучи М. Равновесие в бескоалиционной игре п лиц с выбором момента времен. — С. 55—78.

^ Наумова Н.И. Ограниченная согласованность, порожденная функциями полезности коалиций. — С. 79—99.

^ ПетросянЛ.А., Зенкевич Н.А. Принципы устойчивой кооперации. — С. 100—120.

^ ПетросянЛ.А., Седаков А.А. Многошаговые сетевые игры с полной информацией. — С. 121—138.

^ Тур А.В. Линейно-квадратичные неантагонистические дискретные игры. — С. 139—163.

^ Чуйко Ю.В. Задача маршрутизации с разделяемым трафиком и неполной информацией. — С. 164—176.

^ Галегов А.И., Гарнаев А.Ю. Налоговая игра в дуополии Курно. — С. 177—192.

^ Гарнаев А. Ю., Торицын А.О. Игровая задача справедливого распределения ресурсов при наличии активных помех. — С. 193—208.

^ ГубановД.А., НовиковД.А., Чхартишвили А. Г. Модели репутации и информационного управления в социальных сетях. — С. 209—234. ^ Зенкевич Н.А., Колабутин Н.В., Янг Д.В.К. Стохастическая модель устойчивого совместного предприятия. — С. 235—269.

^ Ивашко А.А. Игра наилучшего выбора двух объектов с полной информацией. — С. 270—286.

^ Искаков М.Б., Павлов П.А. Равновесие в безопасных стратегиях в модели пространственной конкуренции Хотеллинга. — С. 287—318.

^ Коргин Н.А. Эквивалентность и неманипулируемость неанонимных приоритетных механизмов распределения ресурсов. — С. 319—347.

^ Угольницкий Г.А. Оптимизационные и теоретико-игровые модели управления инвестиционно-строительными проектами. — С. 348—365.

^ Реттиева А.Н. Кооперативное регулирующее условие в задаче разделения биоресурсов. — С. 366—384.

^ Шевкопляс Е.В. Уравнение Гамильтона—Якоби—Беллмана в дифференциальных играх со случайной продолжительностью. — С. 385—408.

5.5.1. Качество регулирования

5.5.2. Выбор параметров регулятора

5.5.3. Ручная настройка, основанная на правилах

5.5.4. Методы оптимизации

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

Поэтому для полного описания или тестирования системы с ПИД-регулятором нужен ряд дополнительных показателей качества, о которых речь пойдет ниже.

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

5.5.1. Качество регулирования

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

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

Для той или иной задачи наиболее важными могут быть следующие факторы:

  • форма отклика на внешнее возмущение (время установления, перерегулирование, коэффициент затухания и др.);
  • форма отклика на шумы измерений;
  • форма отклика на сигнал уставки;
  • робастность по отношению к разбросу параметров объекта управления;
  • требования к экономии энергии в управляемой системе;
  • минимум шумов измерений и др.

Для классического ПИД-регулятора параметры, которые являются наилучшими для слежения за уставкой, в общем случае отличаются от параметров, наилучших для ослабления влияния внешних возмущений. Для того, чтобы оба параметра одновременно были оптимальными, необходимо использовать ПИД-регуляторы с двумя степенями свободы (см. раздел «Принцип разомкнутого управления»).

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

Ослабление влияния внешних возмущений

Как было показано в разделе «Запас устойчивости и робастность», обратная связь ослабляет влияние внешних возмущений в раз за исключением тех частот, на которых . Внешние возмущения могут быть приложены к объекту в самых разных его частях, однако, когда конкретное место неизвестно, считают, что возмущение воздействует на вход объекта. В этом случае отклик системы на внешние возмущения определяется передаточной функцией (см. (5.42))

.

(5.109)

Поскольку внешние возмущения обычно лежат в низкочастотной части спектра, где , и, следовательно, , то предыдущее выражение можно упростить:

.

(5.110)

Таким образом, для ослабления влияния внешних возмущений (в частности, влияния нагрузки) можно уменьшить постоянную интегрирований .

Во временной области реакцию на внешние возмущения оценивают по отклику на единичный скачок (см. рис. 5.56).

Ослабление влияния шумов измерений

Передаточная функция от точки приложения шума (рис. 5.35) на выход системы имеет вид (см. (5.42)):

.

(5.111)

Благодаря спаду АЧХ объекта на высоких частотах функция чувствительности стремится к 1 (см. рис. 5.81). Поэтому ослабить влияние шумов измерений с помощью обратной связи невозможно. Однако эти шумы легко устраняются применением фильтров нижних частот, а также правильным экранированием и заземлением [Денисенко, Денисенко].

Робастность к вариации параметров объекта

Замкнутая система остается устойчивой при изменении параметров объекта на величину , если выполняется условие (5.100).

Критерии качества во временной области

Для оценки качества регулирования в замкнутой системе с ПИД-регулятором обычно используют ступенчатое входное воздействие и ряд критериев для описания формы переходного процесса (рис. 5.84):

и момент времени ,
при котором ошибка достигает этого максимума;

(5.112)

;

(5.113)

;

(5.114)

.

(5.115)

Отметим, что в литературе встречаются и другие определения декремента затухания, в частности, как или как коэффициент в показателе степени экспоненты, описывающей огибающую затухающих колебаний;

Рис. 5.84. Критерии качества регулирования во временной области

Рис. 5.85. Критерии качества регулирования в частотной области

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

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

Частотные критерии качества

В частотной области обычно используются следующие критерии, получаемые из графика амплитудно-частотной характеристики замкнутой системы (см. рис. 5.85):

, типовыми значениями являются =1,5…1,6;

(5.114)

Частотные критерии у реальных регуляторов не могут быть однозначно связаны с временными критериями из-за нелинейностей (обычно это нелинейности типа ограничений) и алгоритмов устранения эффекта интегрального насыщения. Однако приближенно можно установить следующие зависимости между критериями в частотной и временной области:

5.5.2. Выбор параметров регулятора

В общей теории автоматического управления структура регулятора выбирается исходя из модели объекта управления. При этом более сложным объектам управления соответствуют более сложные регуляторы. В нашем же случае структура регулятора уже задана — мы рассматриваем ПИД-регулятор, причем эта структура очень простая. Поэтому ПИД-регулятор не всегда может дать хорошее качество регулирования, хотя в подавляющем большинстве приложений в промышленности применяются именно ПИД-регуляторы.

Впервые методику расчета параметров ПИД-регуляторы предложили Зиглер и Никольс в 1942 году [Ziegler]. Эта методика очень проста и дает не очень хорошие результаты. Тем не менее, она до сих пор часто используется на практике, хотя с тех пор появилось множество более точных методов.

После расчета параметров регулятора обычно требуется его ручная подстройка для улучшения качества регулирования. Для этого используется ряд правил, хорошо обоснованных теоретически.

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

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

Настройка параметров регулятора по методу Зиглера и Никольса

Зиглер и Никольс предложили два метода настройки ПИД-регуляторов [Ziegler]. Один из них основан на параметрах отклика объекта на единичный скачок; второй метод основан на частотных характеристиках объекта управления.

Табл. 27. Формулы для расчета коэффициентов регулятора по методу Зиглера-Никольса

Расчет по отклику на скачок

Расчет по частотным параметрам

Регулятор

П

ПИ

ПИД

Примечание. Система обозначений параметров регулятора и формулы соответствует уравнению (5.36).

Для расчета параметров ПИД-регулятора по первому методу Зиглера-Никольса используются всего два параметра: и (см. рис. 5.29 и пояснения к нему в тексте). Формулы для расчета коэффициентов ПИД-регулятора сведены в табл. 27.

Рис. 5.86. Результат настройки ПИД-регулятора по методу Зиглера-Никольса для объекта второго порядка с задержкой: , .

В качестве примера на рис. 5.86 приведен отклик на единичный скачок системы с объектом второго порядка и ПИД-регулятором, настроенным по табл. 27 и переходная характеристика самого объекта управления. Из характеристики объекта получены значения и . По табл.1 для этих значений и можно найти коэффициенты ПИД регулятора: , , . На рис. 5.86 приведен также отклик на единичный скачок той же системы при параметрах , , , полученных путем ручной подстройки. Как видим, метод Зиглера-Никольса дает параметры, далекие от оптимальных. Это объясняется не только упрощенностью самого метода (он использует только 2 параметра для описания объекта), но и тем, что параметры регулятора в этом методе определялись Зиглером и Никольсом исходя из требования к декременту затухания, равному 4, что и дает медленное затухание процесса колебаний.

Метод Зиглера-Никольса никак не учитывает требования к запасу устойчивости системы, что является вторым его недостатком. Судя по медленному затуханию переходного процесса в системе, этот метод дает слишком малый запас устойчивости.

Второй метод Зиглера-Никольса (частотный метод) в качестве исходных данных для расчета использует частоту , на которой сдвиг фаз в разомкнутом контуре достигает 180˚, и модуль коэффициента передачи объекта на этой частоте . О методике определении этих параметров см. раздел «Частотная идентификация в режиме релейного регулирования». Зная параметр , сначала находят период собственных колебаний системы , затем по табл.1 определяют параметры регулятора. Точность настройки регулятора и недостатки обоих методов Зиглера-Никольса одинаковы.

Метод CHR

В отличие от Зиглера и Никольса, которые использовали в качестве критерия качества настройки декремент затухания, равный 4, Chien, Hrones и Reswick (CHR) [Chien] использовали критерий максимальной скорости нарастания при отсутствии перерегулирования или при наличии не более чем 20%-ного перерегулирования. Такой критерий позволяет получить больший запас устойчивости, чем в методе Зиглера-Никольса.

CHR метод дает две разные системы параметров регулятора. Одна из них получена при наблюдении отклика на изменение уставки (табл. 28), вторая — при наблюдении отклика на внешние возмущения (табл. 29). Какую систему параметров выбирать — зависит от того, что важнее для конкретного регулятора: качество регулирования при изменении уставки, или ослабление внешних воздействий. Если же важно и то, и другое, то необходимо использовать регуляторы с двумя степенями свободы (см. раздел «Принцип разомкнутого управления»).

Метод CHR использует аппроксимацию объекта моделью первого порядка с задержкой (5.5).

Табл. 28. Формулы для расчета коэффициентов регулятора по методу CHR, по отклику на изменение уставки

Без перерегулирования

С 20%-ным перерегулированием

Регулятор

П

ПИ

ПИД

Примечание. Система обозначений параметров регулятора и формулы соответствует уравнению (5.36).

Табл. 29. Формулы для расчета коэффициентов регулятора по методу CHR, по отклику на внешние возмущения

Без перерегулирования

С 20%-ным перерегулированием

Регулятор

П

ПИ

ПИД

Примечание. Система обозначений параметров регулятора и формулы соответствует уравнению (5.36).

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

Обратим внимание, что пропорциональный коэффициент в методе CHR меньше, чем в методе Зиглера-Никольса.

5.5.3. Ручная настройка, основанная на правилах

Расчет параметров по формулам не может дать оптимальной настройки регулятора, поскольку аналитически полученные результаты основываются на сильно упрощенных моделях объекта. В частности, в них не учитывается всегда присутствующая нелинейность типа «ограничение» для управляющего воздействия (см. раздел «Интегральное насыщение»). Кроме того, модели используют параметры, идентифицированные с некоторой погрешностью. Поэтому после расчета параметров регулятора желательно сделать его подстройку. Подстройку можно выполнить на основе правил, которые используются для ручной настройки. Эти правила получены из опыта, теоретического анализа и численных экспериментов. Они сводятся к следующему [Astrom]:

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

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

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

Отметим, что применение правил возможно только после предварительной настройки регулятора по формулам. Попытки настроить регулятор без начального приближенного расчета коэффициентов могут быть безуспешными. Сформулированные выше правила справедливы только в окрестности оптимальной настройки регулятора. Вдали от нее эффекты могут быть иными, см. раздел «Классический ПИД-регулятор»

При регулировке тепловых процессов настройка по правилам может занять недопустимо много времени.

5.5.4. Методы оптимизации

Методы оптимизации для нахождения параметров регулятора концептуально очень просты и аналогичны численным методам идентификации параметров объекта (см. раздел «Методы минимизации критериальной функции»). Выбирается критерий минимизации, в качестве которого может быть один из показателей качества или комплексный критерий, составленный из нескольких показателей с разными весовыми коэффициентами. К критерию добавляются ограничения, накладываемые требованиями робастности. Таким путем получается критериальная функция, зависящая от параметров ПИД-регулятора. Далее используются численные методы минимизации критериальной функции с заданными ограничениями, которые и позволяют найти искомые параметры ПИД-регулятора.

Методы, основанные на оптимизации, имеют следующие достоинства:

  • позволяют получить оптимальные значения параметров, не требующие дальнейшей подстройки;
  • не требуют упрощения модели объекта, модель может быть как угодно сложной;
  • позволяют быстро достичь конечного результата (избежать процедуры длительной подстройки параметров).

Однако реализация данного подхода связана с большими проблемами, которые не один десяток лет являются предметов научных исследований. К этим проблемам относится:

  • низкая надежность метода (во многих случаях вычислительный процесс может расходиться и искомые коэффициенты не будут найдены);
  • низкая скорость поиска минимума для овражных функций и функций с несколькими минимумами.

Тем не менее, методы оптимизации являются мощным средством настройки ПИД-регуляторов с помощью специально разработанных для этого компьютерных программ (см. раздел «Программные средства настройки»).

  • Ошибка регулирования давления топлива
  • Ошибка регулирования давления al4
  • Ошибка регистров пользователя атол
  • Ошибка регистрация пользователей на настоящий момент не разрешается
  • Ошибка регистрация не может быть завершена вайбер