Міністерство освіти і науки, молоді та спорту України
Національний університет "Львівська політехніка"
Кафедра КСА
Курсова робота
з навчальної дисципліни :
“Моделювання процесів та елементів систем керування”
Однофазний випрямляч з множенням напруги вдвiчi
Варіант №10б1
Технічне завдання:
1. Навести принципову схему пристрою і описати його роботу.
2. Записати рівняння стану електричного і маґнетного кіл пристрою.
3. Звести отримані рівняння до нормальної форми Коші, записавши їх відносно робочих потокозчеплень, струмів і напруг на конденсаторах.
4. Описати метод Рунге-Кутта четвертого порядку точності з постійним кроком інтегрування.
5. Описати ітераційний метод аналізу усталених режимів :
а) екстраполяційний ε-алґоритм
б) метод середніх значень
в) метод Ньютона.
6. Скласти алгоритм розрахунку перехідних і усталених режимів роботи пристрою.
7. Скласти програму розрахунку перехідних і усталених режимів роботи пристрою алгоритмічною мовою С++.
8. Виконати розрахунок 5 періодів перехідного процесу, навести графічні результати.
9. Виконати розрахунок усталеного режиму, навести графічні результати.
10. За отриманими результатами зробити висновки.
Параметри системи :
Варіант
1, Гн-1
2, Гн-1
r1, Ом
r2, Ом
C1, мФ
C2, мФ
1
150
120
13
11
4
0,4
Варіант
і1, А
i2, A
1, Вб
2, Вб
m2, Гн-1
Um, В
1
0,05
0,9
0,2
0,9
3,0
250
Параметри навантаження:
Варіант
схема а)
схема б)
RН, Ом
RН, Ом
LН, мГн
1
15
13
0,5
Позначення: α1, α2 – обернені індуктивності розсіювання обмоток трансформатора; r1, r2 – активні опори цих обмоток
Однофазний випрямляч з множенням напруги вдвiчi
Принцип роботи пристрою:
Схеми множення напруги застосовують у тих випадках, коли необхiдно мати високу напругу на навантаженнi при заданiй вхiднiй напрузi. Їх доцiльно використовувати при високiй напрузi та малих струмах навантаження.
Подвоєння напруги відбувається таким чином. При від’ємній півхвилі напруги на виході трансформатора відкривається діод і конденсатор заряджається. Зміна полярності вихідної напруги приводить до закривання діоду і відкривання діоду . Тоді починається заряджатися конденсатор , причому, до нього буде прикладена напруга, яка дорівнює сумі вихідної напруги трансформатора та напруги конденсатора . Ця сума наближено буде вдвічі більшою від вихідної напруги. Таким чином, до навантаження буде прикладена подвійна напруга.
Моделювання таких схем виконуємо шляхом введення в рiвняння динамiки пристрою додаткових логiчних змiнних , якi набувають значень 0, 1. Значення 0 вiдповiдає закритому стану діода, а 1 – вiдкритому. Це дозволяє описати єдиною системою рiвнянь всi можливi комбiнацiї вiдкритих i закритих діодів. Таких комбінацій буде дві.
Схеми пристрою:
Принципова електрична схема:
Схема заміщення при закритому діоді Д1 та закритому Д2 (k1= k2=0):
Схема заміщення при відкритому діоді Д1, та закритому Д2 (k1=1, k2=0):
Схема заміщення при відкритому діоді Д2, та закритому Д1 (k1=0, k2=1):
Розглянемо перший випадок, коли обидва діоди закритi.
Схема заміщення при закритому діоді Д1 та закритому Д2 (k1= k2=0):
Розрахункова схема вiдповiдає неробочому ходу трансформатора (=0), тому тут ми будемо мати лише одне диференцiальне рівняння:
де - робоче потокозчеплення; - повне потокозчеплення, напруга, опiр i струм первинної обмотки трансформатора; – величина, обернена індуктивності розсіяння первинної обмотки;
– величина, обернена диференціальній індуктивності трансформатора.
Згідно розрахункової схеми, запишемо рівняння конденсаторів:
Та рівняння котушки:
Струм первинної обмотки можна визначити аналітично з рiвняння стану маґнетного кола
.
де – величина, обернена статичній індуктивності трансформатора
Струм навантаження та другого конденсатора , згідно розрахункової схеми, вiдповiдно рiвнi
Рівняння струму первинної обмотки можна записати в диференціальній формі:
.
Запишемо загальним матричним виразом:
де – вектор змінних стану;
– вектор часових функцій
;
– матриця коефіцієнтів.
Отримані рівняння необхідно доповнити умовами закривання і відкривання діодів. Обидва діоди закриті, тому необхідно контролювати умови їх вiдкривання
,
.
При виконаннi першої умови вiдкривається діод і логічна змінна , при виконанні другої умови – вiдкривається діод і логічна змінна .
Розглянемо рівняння динаміки другої комбінації, коли діод вiдкритий, а – закритий:
При такiй комбiнацiї трансформатор навантажений конденсатором , а конденсатор розряджається на навантаження.
де – величина, обернена індуктивності розсіяння вторинної обмотки трансформатора;
Отримані рівняння необхідно доповнити рівнянням конденсатора
де – величина, обернена індуктивності розсіяння вторинної обмотки трансформатора.
Рiвняння конденсатора :
Та рівняння котушки:
а струм навантаження:
.
Струм первинної обмотки можна визначити аналітично згідно
Запишемо загальним матричним виразом:
– вектор змінних стану;
– вектор часових функцій
; – матриця коефіцієнтів
Отримані рівняння необхідно доповнити умовами закривання і відкривання діодів. Діод вiдкритий, тому необхідно контролювати умову його закривання
Діод закритий, тому необхідно контролювати умову його вiдкривання
.
При виконаннi першої умови діод закривається, а при виконанні другої – вiдкривається діод .
Розглянемо рівняння динаміки третьої комбінації, коли діод закритий, а – вiдкритий :
Тут у вихідному контурі трансформатора добавиться лише ще одна напруга конденсатора , тому рівняння набудуть вигляду:
Рівняння конденсатора , згідно розрахункової схеми, набуде вигляду:
Струм навантаження та конденсатора будуть
Систему нелінійних диференціальних запишемо одним матричним виразом в нормальній формі Коші, де
– вектор змінних стану;
– вектор часових функцій
– матриця коефіцієнтів
;
Отримані рівняння необхідно доповнити умовами закривання і відкривання діодів. Діод вiдкритий, тому необхідно контролювати умову його закривання
Діод закритий, тому необхідно контролювати умову його вiдкривання
.
Всі три описанi комбiнацiї вiдкритих i закритих діодів можна звести до однiєї системи нелінійних диференціальних рiвнянь, записаної в нормальній формі Коші, ввiвши додатковi логiчнi змiннi i , що набувають значень 0, 1.
де
Аналогiчно можна об'єднати рiвняння конденсаторів:
Струми визначаємо за формулами
Запишемо рівняння в нормальній формі Коші, де
– вектор змінних стану;
– вектор часових функцій
.
– матриця коефіцієнтів,
Метод Рунге-Кутта:
Метод Рунге-Кутта об’єднує ціле сімейство методів розв’язування диференціальних рівнянь першого порядку. Найбільш часто використовується метод четвертого порядку, який просто називають “ методом Рунге-Кутта”.
В методі Рунге-Кутта значення функції , як і в методі Ейлера, визначається за формулою :
Якщо розкласти функцію в ряд Тейлора і обмежитись членами до включно, то приріст можна записати у вигляді :
Замість того, щоб обчислювати члени ряду за формулою (2) в методі Рунге-Кутта використовують наступні формули.
Це метод четвертого порядку точності.
Похибка на кожному кроці має порядок .Таким чином метод Рунге-Кутта забезпечує значно вищу точність ніж метод Ейлера, однак вимагає більшого об’єму обчислень ніж метод Ейлера. Це досить часто дозволяє збільшити крок .
Деколи зустрічається інша форма представлення методу Рунге-Кутта 4-го порядку точності.
Екстраполяційний ε-алґоритм:
Екстраполяційні методи цільову функцію використовують у вигляді Інтеґруючи на періодах, породжуємо послідовність значень змінних стану :
(1)
де – період вхідної дії. Для послідовності (1), починаючи з , застосовуємо екстраполяційну формулу :
,(2)
де – початкові умови входження в усталений режим.
В якості функції EXTR доцільно використати – алґоритм, який виконує обчислення границі послідовності з експоненційними складовими. Формула для обчислення наступного значення має вигляд :
(2)
де Результат екстраполяції, що відповідає EXTR в (2) рівний В (3) використовується процедура обертання Самельсона де – k-й елемент n-мірної колонки .
Для систем розмірності значення . Початкові умови усталеного режиму для швидкозатухаючих перехідних реакцій компонент визначаємо інтеґруванням рівнянь стану системи на періодах. Як правило . На жаль, не існує строгого критерію вибору i , тому тут можливий лише еврістичний підхід. Та й визначити які саме компоненти будуть швидкозатухаючими без чисельного експерименту неможливо. Таким чином, основний недолік екстраполяційних методів полягає в необхідності інтеґрування рівнянь динаміки на значному інтервалі часу.
Розрахунок перехідного процесу:
Розрахунок усталеного режиму:
1. Інтеґруємо рівняння , від заданих початкових умов на періодах і визначаємо початкові умови періодичного режиму швидкозатухаючих компонент змінних стану .
2. Маючи на -й ітерації формули початкові умови , інтеґруємо рівняння стану на періодах і породжуємо послідовність :
. (2.30)
3.Згідно
визначаємо уточнене значення початкових умов
(2.31)
4. Перевіряємо умову
збіжності ітераційного процесу. Якщо вона не виконується, то процес повторюємо з п.2, в іншому випадку зупиняємо ітераційний процес.
Текст програми:
//---------------------------------------------------------------------------
#include <vcl.h>
#pragma hdrstop
#include <math.h>
#include <fstream.h>
#include <iostream.h>
#include "Unit1.h"
//---------------------------------------------------------------------------
#pragma package(smart_init)
#pragma resource "*.dfm"
TForm1 *Form1;
using namespace std;
//---------------------------------------------------------------------------
//---------------------------------------------------------------------------
__fastcall TForm1::TForm1(TComponent* Owner)
: TForm(Owner)
{
}
//---------------------------------------------------------------------------
//Початок. Параметри системи-------------------------------------------------
const double Freq=50,T=1/Freq, Um=250;
const double a1=150,a2=120, r1=13,r2=11,C1=4e-3,C2=0.4e-3,Ln=0.5e-3,Rn=13;
const double f1=0.2,f2=0.9,i1=0.05,i2=0.9,m2=3; //ïàðàìåòðè êðèâî¯ íàìàãí
int K_1,K_2;
double I_f (double f);
double I_f (double f)
{
double b1,b2,b3,b4,m1,h;
int i;
h=f2-f1; m1=i1/f1;
if (f<0) i=-1;
else i=1;
f=fabs(f);
if (f<=f1)
return i*(m1*f);
if (f>=f2)
return i*(m2*f+i2-m2*f2);
b1=(2*(f-f1)+h)*(f2-f)*(f2-f);
b2=(2*(f2-f)+h)*(f-f1)*(f-f1);
b3=(f-f1)*(f2-f)*(f2-f);
b4=(f-f2)*(f-f1)*(f-f1);
return i*( (b1*i1+b2*i2)/h/h/h+(b3*m1+b4*m2)/h/h );
}
double a_ii (double f);
double a_ii (double f)
{
double db1,db2,db3,db4,m1,h;
h=f2-f1; m1=i1/f1;
f=fabs(f);
if (f<=f1)
return m1;
if (f>=f2)
return m2;
db1=2*(f-f2)*(3*f-2*f1-f2+h);
db2=-2*(f-f1)*(3*f-f1-2*f2-h);
db3=(f-f2)*(3*f-2*f1-f2);
db4=(f-f1)*(3*f-f1-2*f2);
return (db1*i1+db2*i2)/h/h/h+(db3*m1+db4*m2)/h/h;
}
void D_Control(double *X, double *k4);
void D_Control(double *X, double *k4)
{
//Оголошення зовнішніх змінних
extern int K_1,K_2;
//Аналіз умови відкриття вентилів. Випадок 1
//Вентиль 1
if((-k4[0] - X[2]) > 0)
{
K_1 = 1;
}
//Вентиль 2
if((k4[0] + X[2] - X[3]) > 0)
{
K_2 = 1;
}
//Аналіз умови відкриття-закриття вентилів. Випадок 2
//Закриття вентилів
//Вентиль 1
if(X[1] < 0)
{
K_1 = 0;
X[1] = 0;
}
//Відкриття вентилів
//Вентиль 2
if((X[3] < 0)&&(X[1] < 0))
{
K_2 = 1;
}
//Аналіз умови відкриття-закриття вентилів. Випадок 3
//Закриття вентилів
//Вентиль 2
if(X[1] > 0)
{
K_2 = 0;
X[1] = 0;
}
//Відкриття вентилів
//Вентиль 1
if((X[4] < 0)&&(X[1] > 0))
{
K_1 = 1;
}
}
double dFdt (double *K, double *X, double t);
double dFdt (double *K, double *X, double t)
{
//Розмірність матриці X
const int size = 5;
//Перевизначення вхідних даних
double f=X[0], i2=X[1], Uc1=X[2], Uc2=X[3], in=X[4], i1=a_ii(f)*f-i2;
double al_zn=a1+(K_1+K_2)*a2+a_ii(f);
double g1=a1/al_zn, g2=(K_1+K_2)*a2/al_zn;
double A2[2]={-a2*g1, (K_1+K_2)*a2*(1-g2)};
double U1=Um*sin(2*M_PI*Freq*t);
double dYdt[size]={U1-r1*i1,Uc2*K_2-Uc1-r2*i2, i2,-K_2*i2-Uc2/Rn,Uc2-in*Rn};
double B[size][size]={{g1, g2, 0, 0, 0},
{A2[0], A2[1], 0, 0, 0},
{0, 0, 1/C1, 0, 0},
{0, 0, 0, 1/C2, 0},
{0, 0, 0, 0, 1/Ln}};
for (int i=0; i<size; i++)
{
K[i]=0;
for (int k=0; k<size; k++)
K[i]+=B[i][k]*dYdt[k]; //K=B*dYdt
}
return i1;
}
void RungeKutta (double*X);
void RungeKutta (double*X)
{
//Розмірність системи
const int size = 5;
double k1[size], k2[size], k3[size], k4[size], z[size], h=0.00001;
for (double t=0; t<=4*T; t+=h)
{
dFdt(k1,X,t);
for ( int i=0; i<size; i++)
z[i]=X[i]+k1[i]*h/2;
dFdt(k2,z,t+h/2);
for ( int i=0; i<size; i++)
z[i]=X[i]+k2[i]*h/2;
dFdt(k3,z,t+h/2);
for ( int i=0; i<size; i++)
z[i]=X[i]+k3[i]*h;
dFdt(k4,z,t+h);
for ( int i=0; i<size; i++)
X[i]=X[i]+ h/6*(k1[i]+2*k2[i]+2*k3[i]+k4[i]);
D_Control(X,k4);
}
}
void RungeKutta_Pereh(int S, double *X, TLineSeries *Series1,TLineSeries *Series2,TLineSeries *Series3,
TLineSeries *Series4,TLineSeries *Series5, TLineSeries *Series6);
void RungeKutta_Pereh(int S, double *X, TLineSeries *Series1,TLineSeries *Series2,TLineSeries *Series3,
TLineSeries *Series4,TLineSeries *Series5, TLineSeries *Series6)
{
//Створюємо та відкриваємо файл для запису даних--------------------
FILE *pFile;
if(S == 1)
{
pFile = fopen("KURS_WORK_T_psi_i1_i2_Uc1_Uc2_in_t.dat","w+");
}
else
{
pFile = fopen("KURS_WORK_S_psi_i1_i2_Uc1_Uc2_in_t.dat","w+");
}
//Заголовки даних psi,i1,i2,Uc1,Uc2,in
fprintf(pFile,"%-15.10s %-15.10s %-15.10s %-12.10s %-10.10s %-10.10s %-10.10s\n","psi","i1","i2","Uc1","Uc2","in","t");
double i1_1;
//Розмірність стстеми
const int size = 5;
double k1[size], k2[size], k3[size], k4[size], z[size], h=0.00001;
for (double t=0; t<=5*T; t+=h)
{
i1_1 = dFdt(k1,X,t);
for ( int i=0; i<size; i++)
z[i]=X[i]+k1[i]*h/2;
i1_1 = dFdt(k2,z,t+h/2);
for ( int i=0; i<size; i++)
z[i]=X[i]+k2[i]*h/2;
i1_1 = dFdt(k3,z,t+h/2);
for ( int i=0; i<size; i++)
z[i]=X[i]+k3[i]*h;
i1_1 = dFdt(k4,z,t+h);
for ( int i=0; i<size; i++)
X[i]=X[i]+ h/6*(k1[i]+2*k2[i]+2*k3[i]+k4[i]);
//аналіз умови відкриття-закриття вентилів-------------------------
D_Control(X,k4);
//будуємо графіки psi,i1,i2,Uc1,Uc2,in--------------------
Series1->AddXY(t,X[0],"",clRed);
Series2->AddXY(t,i1_1,"",clRed);
Series3->AddXY(t,X[1],"",clRed);
Series4->AddXY(t,X[2],"",clRed);
Series5->AddXY(t,X[3],"",clRed);
Series6->AddXY(t,X[4],"",clRed);
fprintf(pFile,"%-10.10f %-10.10f %-10.10f %-10.10f %-10.10f %-10.10f %-10.10f\n",X[0],i1_1,X[1],X[2],X[3],X[4],t);
}
fclose(pFile);
}
double* E_algorithm();
double* E_algorithm()
{
//Параметри системи
//m –розмірність системи рівнянь
//n = 2*q + 1 – розмір Sn послідовності
//p – кількість початкових операцій
const int m=5, q=5, n=2*q+1, p=2;
double X[m]={0,0,0,0,0}, e[n+1][n+1][m], eps=0.00001;
double V[m], sum;
bool cond;
int count=0;
//Встановлення нулів для випадку k = -1
for (int j=0; j<=n; j++)
for (int i=0; i<m; i++)
e[j][0][i]=0;
do
{
cond=false;
for (int j=1; j<=p; j++)
RungeKutta (X);
for (int i=0; i<m; i++)
e[0][1][i]=X[i];
for (int j=0; j<=2*q-1; j++)
{
RungeKutta(X);
for (int i=0; i<m; i++)
e[j+1][1][i]=X[i];
if (j==0)
{
for (int i=0; i<=m-1; i++)
{
if(e[1][1][i] == 0) {continue;}
cond=cond || (fabs ((e[1][1][i])-e[0][1][i])/e[1][1][i])*100 > eps;
}
if (cond==false) break;
}
}
if (cond)
{
for (int k=1; k<=n-1; k++)
for (int j=0; j<=n-k-1; j++)
{
for (int i=0;i<=m-1;i++)
V[i]=e[j+1][k][i]-e[j][k][i];
sum=0;
for (int i=0;i<=m-1;i++)
sum=sum+V[i]*V[i];
for (int i=0;i<=m-1;i++)
{
if(sum == 0){continue;}
V[i]=V[i]/sum;
}
for (int i=0; i<=m-1;i++)
e[j][k+1][i]=e[j+1][k-1][i]+V[i];
}
for (int i=0; i<=m-1;i++)
X[i]=e[0][n][i];
count++;
}
}
while (cond);
//Створення динамічного масиву для передачі значень
double *Matrix=(double *)malloc(m*sizeof(double));
//Копіювання значень у динамічний масив
for (int i=0;i<m;i++)
{
Matrix[i] = X[i];
}
return Matrix;
}
void __fastcall TForm1::Button1Click(TObject *Sender)
{
//Вмикаємо кнопку
Button1->Enabled = false;
extern int K_1,K_2;
//початкові умови. Вентилі закриті
K_1 = 0;
K_2 = 0;
//Очищаємо область графіків
Series1->Clear();
Series2->Clear();
Series3->Clear();
Series4->Clear();
Series5->Clear();
Series6->Clear();
int Switch = 0;
//матриця початкових умов, що задовольняють усталеному процесу
double *X;
//Розрахунок усталеного режиму
X = E_algorithm();
//моделювання усталеного режиму
RungeKutta_Pereh(Switch,X,Series1,Series2,Series3,Series4,Series5,Series6);
//Вивільнення памяті
free(X);
//Вмикаємо кнопку
Button1->Enabled = true;
}
void __fastcall TForm1::Button2Click(TObject *Sender)
{
//Розмірність системи
const int size = 5;
//Вмикаємо кнопку
Button2->Enabled = false;
extern int K_1,K_2;
//Початкові умови. Вентилі закриті
K_1 = 0;
K_2 = 0;
Series1->Clear();
Series2->Clear();
Series3->Clear();
Series4->Clear();
Series5->Clear();
Series6->Clear();
int Switch = 1;
double X[size] = {0,0,0,0,0};
//Моделювання перехідного процесу
RungeKutta_Pereh(Switch,X,Series1,Series2,Series3,Series4,Series5,Series6);
//Вмикаємо кнопку
Button2->Enabled = true;
}
//---------------------------------------------------------------------------
Результати роботи. Перехідний процес:
/
Усталений режим:
/
Висновки:
Виконуючи курсову роботу, я закріпила знання і навики аналізу перехідних і усталених режимів роботи чисельними методами нелінійних електромаґнетних пристроїв систем автоматичного керування.
Для моделювання роботи пристрою були використані метод Рунге-Кутта 4 порядку з постійним кроком, а також екстраполяційний ε – алгоритм. Вибір останнього обумовлений тим, що в дану систему входять напівпровідникові вентилі(діоди). Це призводить до зміни системи рівнянь, які описують пристрій і обчислення матриці монодромії згідно моделі чутливостей до початкових умов стає ускладненим. Цього недоліку позбавлений ε – алгоритм.