49582

РОЗРОБКА МАТЕМАТИЧНИХ МОДЕЛЕЙ ЕЛЕКТРОННИХ СХЕМ У РІЗНИХ РЕЖИМАХ ЇХ РОБОТИ

Курсовая

Коммуникация, связь, радиоэлектроника и цифровые приборы

До основних якісних показників і параметрів підсилювача відносяться коефіцієнт передачі коефіцієнт підсилення Кр вхідний і вихідний опори Zвх Zвих динамічний діапазон коефіцієнт нелінійних викрівлень коефіцієнт шуму...

Украинкский

2014-01-03

1.71 MB

1 чел.

Міністерство освіти України

Національний технічний університет України

«Київський політехнічний інститут»

Інститут телекомунікаційних систем

 

Теорія електричних ланцюгів

Методичні вказівки до роботи

«РОЗРОБКА МАТЕМАТИЧНИХ МОДЕЛЕЙ ЕЛЕКТРОННИХ СХЕМ У РІЗНИХ РЕЖИМАХ ЇХ РОБОТИ»

Б.М.Шелковніков

Розглянені і ухвалені

на засіданні інститута

телекомунікаційних мереж і систем

Протокол №_______________________

      від   ________________________

Київ - 2008

УДК  621.395.001

Теорія електричних ланцюгів

Методичні вказівки до роботи

«РОЗРОБКА МАТЕМАТИЧНИХ МОДЕЛЕЙ ЕЛЕКТРОННИХ СХЕМ

У РІЗНИХ РЕЖИМАХ ЇХ РОБОТИ»

Б.М.Шелковніков

-НТУУ «КПІ», 2008р.

Методичні вказівки вміщують три розділи:

    Математичні моделі різних режимів роботи електронних схем і

    методи та алгоритми розрахунку різних режимів роботи електронних схем.

    Методи и алгоритми аналіза чутливості електронних схем.

    Методи и алгоритми оптимізації електронних схем.

Методичні вказівки до курсової роботи призначені для виконання курсової роботи студентами інститута телекомунікаційних систем  НТУУ «КПІ», з вищезгаданої дисципліни, а також для самостійної роботи при вивченні курса.

Бібліографія назв.

Рецензенти:  

     

ВСТУП

Використання персональних електронних обчислювальних машин (ПЕОМ) у всіх областях людської діяльності - характерна риса науково-технічної революції. ПЕОМ, особливо  високопродуктивні, сприяють прискоренню прогресу у радіоелектронній промисловості. Використання ПЕОМ передбачає розробку відповідного спеціалізованого математичного (методи, алгоритми) і програмного забезпечення.

Мета курсу викладеного в методичних вказівках – допомогти у вивченні електронних схем як об'єктів дослідження і проектування, отримання навичок у  формулюванні задач дослідження і проектування, оволодіння методами і алгоритмами рішення задач дослідження у проектуванні електронних схем, навичками реалізації задач у вигляді програмного забезпечення  на ПЕОМ.  Виклад  курсу базується  на  знаннях студентами курсів з математики, фізики, теоретичних основ електротехніки, напівпровідникових приладів, електронних ланцюгів безперервної та імпульсної дії.

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

ЗАДАЧІ АНАЛІЗУ ЕЛЕКТРОННИХ  СХЕМ

Транзисторний підсилювач (ТРП), представлений електричною принциповою схемою (на мал. 1), в залежності від характера вхідного сигналу може працювати у різних режимах. За відсутності вхідного сигналу (або постійному сигналі) підсилювач знаходиться у статичному стані (режим постійного струму). При малому вхідному сигналі, що швидко змінюється, припустимо вважати, що транзистор має тільки лінійні динамічні властивості, і підсилювач працює в режимі лінійного підсилення. При великому вхідному сигналі, що швидко змінюється, транзистор виявдяє нелінійні динамічні властивості, підсилювач функціонує в динамічному нелінійному режимі. У залежності від форми вхідного сигналу (гармонічний, імпульсний) функціонування підсилювача може розглядатися у тимчасовій або частотній областях.

Мал 1

Кожний режим роботи підсилювача можна представити відповідним еквівалентним ланцюгом (схемою) і математичною моделлю і оцінити багатьма якісними показниками (характеристиками, схемними функціями) і параметрами. Якісні показники визначаються на основі математичної моделі і перевіряються експериментально. Весь перелік якісних показників характеризує властивості і функціональні можливості підсилювача загалом.  До основних якісних показників і параметрів підсилювача відносяться коефіцієнт передачі (коефіцієнт підсилення) Кр, вхідний і вихідний опори Zвх, Zвих,  динамічний діапазон, коефіцієнт нелінійних викрівлень, коефіцієнт шуму. Щоб знайти ці якісні показники необхідно проаналізувати підсилювач у статичному режимі, у динамічному режимі у часовій і частотній областях при великому і малому вхідних сигналах.

Еквівалентна схема ТРП для кожного режиму має свою множність елементів (компонентів) і свою структуру (тобто специфічне для режиму зєднання елементів). Так, наприклад, режим малого вхідного сигналу представляється лінійною еквівалентною схемою – з’єднанням лінійних елементів, статичний режим - нелінійною еквівалентною схемою на постійному струмі і т.д.

Слід зазначити, що згадані режими характеризують роботу більшості електронних схем приймально-підсилювальних пристроїв і тому розв’язання задач розрахунку схем в цих режимах має загальне значення.               

Множина якісних показників, що визначаються у відповідному режимі  і представляють задачі аналізу, залежить від множності елементів і їх параметрів - рэ, від структури їх з'єднання - Sp, типу вхідного сигналу (постійний, частотний, часовий):

Кр=F(Sp, рэ, Up)                       (1)

де p - відповідний режим.

Динамічні якісні показники завжди залежать від вихідного статичного режиму, що можна відобразити залежністю:

рэ=( Кст)

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

Велике значення при проектуванні електронних схем має розв’язання задач розрахунку чутливості якісних показників по параметрах елементів- Sр, що дозволяє визначити допуски на параметри, і задач оптимізації, тобто пошуку множини оптимальних параметрів ропт, що забезпечують необхідне відхилення множини якісних показників від заданих у технічному завданні.

 Серед згаданих режимів найбільш загальним є динамічний режим при дії великого сигналу, що змінювався у часі. Інші режими – часткові випадки від цього режиму.

Динамічний нелінійний режим (часова і частотна область) – при великому вхідному сигналі.

 Статичний режим спостерігається, коли зовнішній вплив постійний у часі.

 Динамічний лінійний режим (часова і частотна область) - при малому вхідному сигналі.

МАТЕМАТИЧНА МОДЕЛЬ ЕЛЕКТРОННОЇ СХЕМИ В ДИНАМІЧНОМУ РЕЖИМІ ПРИ ВЕЛИКОМУ СИГНАЛІ.

ЧАСОВА ОБЛАСТЬ

При великому вхідному сигналі, що швидко змінюється, в ТРП транзистор виявляє нелінійні і динамічні властивості (мал. 2), які можуть бути представлені еквівалентною схемою (на мал. 2 виділена штрихами) і математичною моделлю Еберса-Молла [1,2].

Мал. 2

У моделі Еберса-Молла властивості елементів виражаються наступними співвідношеннями:

Iк=Iко•(e(К1•Uкб)-1)=К(Uкб),  К1=1/(mкТ),                         (2)

Iэ=Iэо•(e(К2•Uбэ)-1)=Э(Uбэ),  К2=1/(mэТ)

Jдк=N•Iэ,

Jдэ=I•Iк,

URб=iRб•Rб

Iск=Ск(Uкб)•dUкб/dt,

Ск(Uкб)=Cкб+Скд=Сокб/(1-Uкб/К)0.5+Сокд•Iк

Сокд=К1/(2••Fi)

Iсэ=Сэ(Uбэ)•dUбэ/dt,

Сэ(Uбэ)=Cэб+Сэд=Сокб/(1-Uбэ/К)0.5окд•Iэ

Соэд=К2/(2••Fn)

де Uкб, Uбэ - напруга колектор-база, база-еміттер відповідно;

К1, К2 - температурний потенціал;

т - контактна різниця потенціалів;

Iко, Iэо – струми насичення колекторного та еміттерного переходів;

mк, mэ – коефіцієнти, що відображують технологію виготовлення транзисторів;

N, I - коефіцієнти підсилення за струмом при нормальному та інверсному режимах;

iструм через резистор бази;

Rб- опір бази;

Сокб, Соэб - бар'єрні ємності при нульовому зміщенні;

Fn Fi - гранична частота транзистора при нормальному та інверсному включеннях.

Властивості інших елементів еквівалентної схеми ТРП (мал. 2) для динамічного режиму виражаються співвідношеннями:

UR1=iR1•R1, UR2=iR2•R2, URб=iRб•Rб, UR3=iR3•R3, UR4=iR4•R4, UR5=iR5•R5,

UR6=iR6•R6,  UR7=iR7•R7

У загальному вигляді можна для елементів схеми можна записати:

Резистори:       URi=iRi•Ri ,                           i - номер резистора                         (3)

Ємності:           iCj=Cj• dUcj/dt ,                  j - номер ємності                  

Джерела постійного струму:  J==const

Вхідний струм:    J~=(t)=Jм•sin(•t+J)- функція часу.

Елементи схеми (або гілки) зєднуючись у вузлах утворюють контура. До струмів у вузлах (перерізах) схеми і напруг у контурах можна застосувати, відповідно, перший і другий закони Кирхгофа:

  1.  Алгебраїчна сума струмів i у будь-якому вузлі (у замкненому перерізі) електричної схеми дорівнює нулю (витікаючий з вузла струм береться зі знаком "+", втікаючий у вузол  струм береться зі знаком "-" ).
  2.  Алгебраїчна сума напруг u гілок у будь-якому контурі електричної схеми дорівнює нулю.

Рівняння з’єднань, складені за законами Кірхгофа, визначаються тільки схемами з’єднань гілок, тобто геометричною структурою ланцюга, і не залежать від виду і характеристик елементів, тобто фізичного змісту гілок. Тому при складанні рівнянь з’єднань зручно відволікатись від виду і характеристик гілок ланцюга і замінювати їх лініями, що з'єднують вузли, із збереженням числа гілок і вузлів. У результаті отримують так званий лінійний граф (топологічний граф), який являє собою сукупність або систему вузлів (вершин), що зображаються точками, і гілок (ребер), що зображаються відрізками ліній, які з’єднують будь-яку пару вузлів. Таким чином, елементами графа є вузол і гілка  (мал. 3).

Мал. 3

Обєднані множини рівнянь гілок (компонентних рівняннь (2), (3) і топологічних рівнянь (4)) складають математичну модель схеми (ММС) для динамічного режиму при великому сигналі. Якщо схема має l гілок, то кількість рівнянь і кількість змінних ММС дорівнює l•2  при виборі незалежних перерізів і контурів. Для нашої схеми при вказаних стрілками напрямах струмів рівняння (4) має вигляд:

Вузол  1 iR1 + iС1 - J~ =0

Вузол 2 -iС1 +iR3 +iR2 +iRб =0

Вузол  3 -iRб +iСк +iСэ - iК - iЭ - iДК - iДЭ =0                                   (5)

Вузол  4 -iR5 - iС2 - iСэ + iЭ + iДЭ =0

Вузол  5 iR4  + iС3 - iСк + iК + iДК =0

Вузол  6 -iR4  - iR2 - iR7 + J= =0

Вузол  7 -iC3  + iC4 - iR6 =0

Крім струмів і напруг гілок, введемо у розгляд нові змінні - потенціали вузлів i відносно базисного вузла (0=0). Як базисний, зручно взяти вузол, спільний для входу і виходу схеми. Тоді згідно з другим законом Кирхгофа, напруги всіх гілок u і вузлові потенціали  i зв’язуються співвідношеннями:

uR1=1-0, uC1=1-2, uR2=2-6, uR3=2-0, uС3=5-7 uRб=2-3, uR4=5-6, uR6=7-0, uR5=0-4, uC4=7-0 uС2=0-4, uCк=3-5, uCэ=3-4, uIк=5-3, uIэ=4-3 uJдк=5-3, uJдэ=4-3                                                        (6)

Множини рівнянь (5) і (6) можна записати у матричній формі в загальному вигляді

                           | A | • | i |=0                              (7)

| u |=| At | • |  |                                      (8)

де | i |=|iR1 , iС1, iR2 ……, J~, J= |t- вектор струмів усіх гілок схеми;

   | u |=|uR1 , uС1, uR2 ……, uJдк , uJдэ |t -вектор напруг усіх гілок;

   |  |=|1,2,3,4,….q |t - вектор вузлових потенціалів;

q – кількість вузлів,  t- знак транспонування.

Матриця |A|, що має назву матриця інциденцій вузол-гілка, для схеми представлена на мал. 3 і характеризує її структурні властивості. Матриці |A| і співвідношенням (7)-(8) відповідає топологічний (направлений) граф схеми, побудований на множині змінних схеми i, u та i. Граф є геометричним представленням структури схеми. На графі виділені вузли 1,2,3,4,5,6,7. Вибір напряму струмів у гілках графа визначає систему незалежних струмів і напруг в ММС. Виразимо рівняння (5) використовуючи рівняння (2),(3) і (6). У результаті отримаємо систему рівнянь (9): 

Узел   1 (1-0)/R1 +С1•d(1-0)/dt - J~ =0

Узел  2 -С1•d(1-0)/dt +(2-0)/R3  +(2-6)/R2 +(2-3)/Rб =0

Узел  3 -(2-3)/Rб +Ск(3-5)•d(3-5)/dt +Сэ(3-4)•d(3-

- 4)/dt- К(5-3)-Э(4-3) -NIэ(4-3) -IIК(5-3)=0

Узел  4 -(0-4)/R5-С2•d(0-4)/dt-Сэ(3-4)•d(3-4)/dt+

+Э(4-3)+ +IIК(5-3)=0                    (9)

Узел  5 (5-6)/R4 +С3•d(5-7)/dt +Ск(3-5)•d(3-

-5)/dt+К(5-3)+N(4-3)=0

Узел  6 -(5-6)/R4 -(2-6)/R2-(6-0)/R7+ J= =0

Узел  7 -С3•d(5-7)/dt +С4•d(7-0)/dt -(7-0)/R6 =0

Ці рівняння, що звуться вузловими, складені методом, подібним методу вузлових потенціалів для лінійних ланцюгів. Система (9) - це система алгебро-диференціальних нелінійних рівнянь відносно змінних 1,2,3,4,5,6,7,J=,J~. Вона складається з трьох груп рівнянь: лінійних алгебраїчних (7), лінійних диференціальних (1,2,6), нелінійних диференціальних (інші рівняння). Відповідно, змінні діляться на лінійні Xл=J= та J~, лінійні диференціальні Xлд=|1,2,6,7|t, нелінійні диференціальні Xнд=|3,4,5|t.

З урахуванням викладеного система (9) може бути записана у скороченому вигляді:

л=( Xл, Xлд, Xнд)=0;   лд=( Xл, X'лд, Xнд)=0

нд=( Xл, Xлд, X'лд, Xнд, X'нд)=0              (10)

де -л - лінійний оператор;  н - нелінійний оператор.

X'нд,  X'лд – похідні змінних за часом.

 

Множину гілок у схемі (2)-(3) відповідно властивостям їх рівнянь, можна розділити на характерні підмножини гілок:                                      (11)

- джерел струму J=, J~;

- лінійних резисторів R і проводимоcтей G:  UR=iR•R  або  iR=UR•G , де G=1/R

- нелінійних резисторів      Iн=(UR);         

- залежних джерел струму     Iд=•Iн ;

- лінійних ємностей      iСЛ= Сл•d(Ucл)/dt ;

- нелінійних ємностей      iСН= Сн( Ucн)•d(Ucн)/dt .

 

Цьому розподілу відповідає розподіл топологічної матриці |A| на субматриці і запис топологічних рівнянь (7)-(8) у формі:

Гілки

Вузли

J=

J~

R1

R2

R3

RБ

R4

R5

R6

R7

Iк

Iэ

Iдк

Iдэ

С1

С2

С3

С4

Ск

Сэ

1

0

-1

1

0

0

0

0

0

0

0

0

0

0

0

1

0

0

0

0

0

2

0

0

0

1

1

1

0

0

0

0

0

0

0

0

-1

0

0

0

0

0

3

0

0

0

0

0

-1

0

0

0

0

-1

-1

-1

-1

0

0

0

0

1

1

4

0

0

0

0

0

0

0

-1

0

0

0

1

0

1

0

-1

0

0

0

-1

5

0

0

0

0

0

0

1

0

0

0

1

0

1

0

0

0

1

0

-1

0

6

1

0

0

-1

0

0

-1

0

0

-1

0

0

0

0

0

0

0

0

0

0

7

0

0

0

0

0

0

0

0

1

0

0

0

0

0

0

0

-1

1

0

0

АЕ

АR

AH

AД

АСЛ

АСН

JЕ

АЕ

АR

AH

AД

АСЛ

АСН

JR

IН

= 0

(12)

IД

iСЛ

uЕ

АtЕ

1

uR

АtR

2

uН

=

AtH

3

(13)

uД

AtД

uСЛ

АtСЛ

uСН

АtСН

q

Підставимо рівняння (13) у співвідношення (11) і результат у (12). Після перетворень отримаємо матричне рівняння:           

                                                                                                                                 (14)

Е|•|iЕ |+|АR|•|G|•|АtR|•||+|АH|•(|АtH|•||)+|АД|•Н(|АtД|•||)+

+|АСЛ|•|СЛ|•d(|Аt|•||)/dt+|АСН|•|СН(|Аt|•||)|•d(|Аt|•||)/dt = 0

 Рівняння (14) - це записане у більш загальній формі (з урахуванням топологічних субматриць) рівняння (9). Підстановка субматриць і рівнянь гілок на основі (11), (2), (3) і подальше перетворення дадуть у кінцевому результаті (9). Рівняння (14) як і (9) можна представити у формі (10).

Отже, рівняння (9), (6), (2), (3) складають математичну модель ТРП у динамічному режимі, а співвідношення (14), (13), (11) - математичну модель для динамічного режиму класу електронних схем, що представлений на основі множини компонентів (гілок) виду (11). Назвемо цю модель ММС-ДР1 (Математична модель схеми - динамічний режим 1).

Розглянемо ще один вид математичної моделі схеми ТРП (мал. 4):

Мал. 4

Вибір дерева на топологічному графі схеми визначає не тільки системи лінійно-незалежних рівнянь, складених по законах Кирхгофа, але й, у кінцевому результаті, вид і властивості математичної моделі схеми.

Найбільш загальні топологічні властивості електронних схем представляються законами Кирхгофа у формі:

i| =0,                                                                 (15)

u|=0,                                                            (16)

де i|, |Рu| - матриці перерізів і контурів.

Насправді,  під |П| та |Р|  надалі маються на увазі матриці головних перерізів і контурів. Якщо узагальнені вузли, що утворюються перерізами графа, співпадають з вершинами (вузлами) графа, то матриця |П| співпадає з |А|.  Це випадок побудови канонічних перетину і дерева. Побудова дерева (мається на увазі фундаментальне дерево) розбиває гілки графа на гілки дерева (ребра) і гілки, що не увійшли в дерево, звані хордами (зв'язками). Рівняння Кирхгофа при цьому набувають вигляду:

iT

Пi

=

1  |

=0  

или

iT=-iX                  (17)

iX

 

uT

Рu

=

 | 1

=0  

или

uХ=uТ                   (18)

uX

У випадку збігу фундаментального дерева з деревом графа виконується співвідношення: =-t, =-t             (19) 

З урахуванням цього співвідношення (17)-(18) запишуться:

iT= -iX,       uХ=tuТ                                         (20)  (21)                                                                              

На основі рівнянь (20)-(21) може бути складена ММС у формі нормальних звичайних диференціальних рівнянь (у вигляді рівнянь змінних стану [1,2,4,5].

Виберемо дерево графа схеми так, щоб у гілки дерева увійшли джерела напруги, ємності і потрібна кількість резисторів, а в хорди – джерела струму, індуктивності (якщо є) і резистори, що залишилися.  Це завжди можна зробити, якщо гілки схеми не утворюють топологічних вироджень (наприклад, контурів з ємносних гілок і джерел ЕРС або перерізів, утворених індуктивними гілками і джерелами струму). Інакше необхідно усунення виродження [5]. На мал. 5 зображено дерево і перерізи на графі схеми ТРП:

Мал. 5

Відповідна матриця перерізів П для схеми:

Перерізи

гілки

E

Uвх

С1

С2

С3

С4

Сэ

Ск

R3

R1

R2

RБ

R4

R5 

R6

R7

Iд 

Iк

Iдэ

Iдк

7

1

1

1

1

1

-1

2

1

-1

5

1

1

-1

1

-1

|П|=

8

1

-1

9

1

1

1

4

1

-1

1

1

-1

1

6

1

1

1

1

-1

3

1

1

-1

-1

Ємнісний контур (С2СэСкС3С4) розірваний включенням невеликого R7. Матриця |П|  розбивається на ряд характерних субматриць

            

1

ERX

EI

CRX

CI

CHIH

CHIД

RBRX

RBI

RX=

ERX

RX=

EI

CRX

CI

RBRX

RBI

Підставимо компонентні рівняння (11) у (20)-(21), у результаті отримаємо:

С•d(Uc)/dt=CRXGХuRX+CIH+IH(uH)+CIДIH                         (22)

iEB=EBRX•GХuRX+EBIIX    

           (23)

iRB=RBRX•GХ•uRX+RBI•IX 

uRX

tRX

EB

                      EB=|E, uBX|t

                      uRB=uR3                          (24)

                      IX=|IH, IД|t

=

uC

uIX

tIX

uRB

    uRB=RBiRB

У цих рівняннях:

C1

1/R1

C2

1/R2

C=

Cэ(uЭБ)

, GX=1/RX=

1/RБ

Cэ(uКБ)

СН

1/R7

RB=R3,

IH=

IЭ(uЭБ)

IК(uКБ)

IД=•IH=

N

I

IЭ

IК

Підстановка рівняння (24) в (22) призводить останнє до вигляду (25):

-1

EB

EB

С•d(Uc)/dt=

С

(CRX•GХtRX

uC

+CIH•IH•(tIX

uC

)+CIД•IH )

uRB

uRB

Співвідношення (25), (23) і (24) - ММС у формі рівнянь змінних стану, назвемо її ММС-ДР2. (Математична модель схеми - динамічний режим 2).

У скороченому вигляді ці співвідношення запишуться:

Xд=д( Xд, Xл, ЕВ),      Xд=| Xлд, Xнд|t,    Xл=л( Xд, Xл, ЕВ)              (26)

Друге рівняння можна подати і у формі:                 |А|•X=д( Xд, ЕВ)

Якщо підставити у (23)-(25) значення топологічних субматриць IJ і параметрів гілок C1,C2,...,R1,R2,..,  і т.п., то отримаємо математичну модель ТРП для динамічного режиму при дії великого сигналу.

Порівняння ММС-ДР1 і ММС-ДР2 показує, що передусім вони відрізняються видом математичних рівнянь (у першому випадку - неявна форма алгебро-диференціальних нелінійних рівнянь, у другому випадку похідні змінних диференціальних рівнянь виражені явно) і кількістю незалежних змінних (у першому випадку - це , у другому - |uС, uRX, uIX |t. У разі виникнення необхідності отримання напруг і струмів усіх гілок у ММС-ДР1 виникає потреба працювати з (2•l+q) рівняннями (14), (13) і (11), у ММС-ДР2 з n 2•l рівняннями.

МАТЕМАТИЧНА МОДЕЛЬ ЕЛЕКТРОННОЇ СХЕМИ У СТАТИЧНОМУ РЕЖИМІ

Схема знаходиться у  статичному режимі, якщо на неї впливають постійні у часі сигнали, тобто при t=to (або рівному нулю)

uвх(to)=Е=const.

При цьому струми в ємностях (напруги на індуктивностях) рівні нулю, що відповідає duC /dt =0 та  diL /dt=0 або відсутності змін струмів і напруг у схемі. Підставляючи ці умови у співвідношення (9), (14), (25) отримаємо відповідні математичні моделі для статичного режиму – ММС ТРП: ММС-Cтl, ММС-Ст2.

Математична модель ТРП для статичного режиму буде мати вигляд (9) без членів з похідними і при uвх.

Рівняння (6) залишається без змін і дозволяє визначити напруги на всіх гілках (включаючи ємнісні) у статиці після знаходження з (27)  і підстановки у (6). Необхідні струми гілок, як і в динаміці, можна знайти з рівнянь (2), (3).

Модель ММС-Ст1 на основі ММС-ДР1 (див. співвідношення (14), (13), (11)) запишеться як

Е|•|iЕ|+|АR|•|G|•|АtR|•||+|АH|•(|АtH|•||)+|АД|•Н(|АtД|•||) = 0    (27)

|u|=|А|t•||,  |u|=|uЕ , uR, uH, uД, uСЛ, uСН |t         (28)

 

R1

iRi= G•u Ri

G=1/R       

R=

RБ

R7

Iн=Н(uН),   Iн=| Iк , Iэ |t,

Iэ=Iн,   Iд=| Iдк , Iдэ |t,  =|n i|         (29)

      

Скорочено рівняння (27) представляються у формі (див. (10))

л( Xл, Xн)=0              (30)

н( Xл, Xн)=0

де л -  лінійний оператор; н - нелінійний оператор.

Xл,  Xн - незалежні змінні, відповідно, лінійних і нелінійних алгебраїчних рівнянь.

Із співвідношень (25), (23), (24) з урахуванням iС=0 и uвх(to)=Е отримаємо математичну модель електронних схем ММС-Ст2

EB

EB

CRX•GХtRX

uC

+CIH•IH•(tIX

uC

)+CIД•IH =0

uRB

uRB

iЕВ= ЕВRX•GХ•uRX + ЕВI•IХ  IХ=| Iн , Iд |t,                   (31)

iRВ= ЕВRX•GХ•uRX + RВI•IХ  uRX=RВiRВ,

uRX

tRX

EB

=

uC

uIX

tIX

uRB

Скорочено ММС-Ст2 має вигляд (30).

Якщо лінійні рівняння розглядати як частковий випадок нелінійних рівнянь, то усі ММС-Ст і рівняння (30) можна подати як одне операторне рівняння:

( X)=0                            (32)

МЕТОДИ І АЛГОРИТМИ  РОЗРАХУНКУ СТАТИЧНОГО РЕЖИМУ РОБОТИ

Взагалі існує два основних підходи до розвязання задачі розрахунку статичного режиму.

Перший заснований на представленні статичного режиму, до якого прямують при t перехідні процеси у схемі при підключенні до неї джерел живлення і вхідного джерела (його постійної складової). При цьому використовуються динамічна математична модель схеми і методи чисельного інтегрування для її розвязання. Другий підхід заснований на розв’язанні алгебро-трансцендентних нелінійних рівнянь із застосуванням ітераційних, проекційних методів, методів спуску і продовження рішення за параметром, комбінованих методів [1,2,4,5]. 

Найбільше поширення при машинному проектуванні електронних схем знайшов метод Ньютона і його модифікації. Нехай задана система нелінійних алгебраїчних і трансцендентних рівнянь вигляду (32) і відомо, що в деякій області G змінних  (X1, Х2,..., Хn) існує єдиний розв’язок (X*1, Х*2,..., Х*n). Метод Ньютона полягає у тому, що за початковим наближенням змінних (X0102,...,Х0n) знаходиться наступне наближення за формулами:

Xi1=Хi0 -|W(Хi0)|-1( Хi0)  i=1,…n

або

W(Хi0)Хi0= -(Хi0),  Хi1= Хi0+Хi0,

де i0) - значення лівої частини системи (32) при Хi0 називається вектором невязок,  

Wi0)=di0)/dХi0матриця Якобі (якобіан) системи (32),

 Хi0- вектор поправок.

Використовуючи знайдені величини розраховуємо

W(Хi1)Х11= -(Хi1),  i=1,…,n

 Хi2= Хi1+Хi1,     і т.д.

Якщо знайдене k-е наближення, то (k+1)-e наближення обчислюється за формулою

W(Хik)Х1k= -ik),     (33)

 Хik+1= Хik+Хik

Якщо Lim(Хik)k для i=1,…,n тобто Хik (похибка), то кажуть, що метод Ньютона збігається до розвязку.

Як видно з (33) на кожній ітерації процесу наближення до розв’язку потрібно обчислювати значення вектора нев’язок ik), Якобіана W=dik)/dХik , розв’язувати систему лінійних алгебраїчних рівнянь відносно вектора поправок Хik і знаходтити наступне наближення Хik+1 через Хik   і  Хik  за формулою додавання векторів.

Наближений розвязок Хik+1= Хi* бажано отримати з наперед заданою точністю. На практиці досягнуту у процесі ітерацій точність оцінюють за нормою вектора поправок Хik  або за нормою вектора невязок[ik)]. Зрозуміло, що при Хik+1Хi* маємо [Хik]0 и [ik)]0. Звідси випливає, що обчислення потрібно припиняти, якщо [Хik ] < або [ik)] < . Під нормою вектора Хik  або ik) може розумітися або Евклідова норма е - норма 

n               

                                         [Хik] =(  (Хik)2)0.5

       i=1             

або S - норма

       n               

                                         [Хik] =    |Хik|

      i=1             

або рівномірна норма ( m - норма)

                            

                                         [Хik] = max |Хik|

       1 i   n             

Швидкість збіжності методу Ньютона квадратична

Хik+1 k( Хik)2  

Якщо помилка Хik мала, наприклад, Хi << 1, то подальша помилка буде зменшуватись до збільшеного у k - раз квадрата попередньої помилки. Після  кожної ітерації спостерігається подвоєння кількості правильних десяткових знаків у результаті. Для збіжності процесу Ньютона до розвязку Х* необхідно, щоб:

а) початкове наближення Х0 було близько задане до коренів Х*;

б) вектор-функція (Х) має бути визначена і безперервна разом зі своїми частковими похідними першого і другого порядку у деякій області;

в) матриця Якобі W(Х) повинна мати обернену обмежену матрицю;

г) матриця других часткових похідних функції (Х) також має бути обмежена. Ці умови математично складні для апріорного визначення факту збіжності і швидкості збіжності, тому не наводиться їх суворе математичне формулювання, але дається пояснення на конкретних прикладах як вони впливають на процес збіжності і результат розвязку.     

  

ПРИКЛАД МАТЕМАТИЧНОЇ МОДЕЛІ ДЛЯ СТАТИЧНОГО РЕЖИМУ РОБОТИ ЕЛЕКТРОННОЇ СХЕМИ

Розглянемо найпростішу математичну модель ланцюга (мал. 6) з діодом для статичного режиму

        Мал. 6

-E+Uд+RIо•(eUд / т-1)=(Uд)=0       (34)

де      Iо,  Т - параметри моделі діода для статичного режиму;

E, R  - параметри ланцюга.

 Якобіан рівняння (34) буде дорівнювати:

W()=d()/d =1+(1/Т) •R•Iо•eUд / т         (35)

Ітераційна формула Ньютона (див. (33)) з урахуванням (34) і (35) набуде вигляду

                              Uкд/Т                                              Uкд/Т

      [ 1+(1/Т) •R•Iо•e              ] •Ukд = E-Uкд -R•Iо•(e          -1)               (36)

Uк+1д = Uкд+Uкд

На мал. 7 зображена геометрична інтерпретація процесу розвязання за ітераційною формулою (36).

Мал. 7

Процес розв’язання розпочинається з початкового наближення U0Д   і закінчується у близькому околі кореня U*Д. Видно, що вибір початкового наближення U0Д праворуч від кореня призводить до закінчення ітераційного процесу у пошуку розв’язку за 3-4 ітерації. Нахил дотичної в точці, наприклад, [U0Д, (U0Д)] визначаєтся W(U0Д), а приріст між ітераціями U0Д – значенням Якобіана і функції у попередній точці, тобто

U0Д = W-1(U0Д) • (U0Д)

У проміжку між ітераціями функція (U0Д) замінюється прямою лінією, дотичною у початковій точці. Тому метод Ньютона називають ще методом дотичної або методом лінеаризації. Якщо початкове наближення вибрати зліва від кореня, то неважко бачити, що через велике значення оберненої похідної W-1(U0Д) вже перший приріст U0Д великий і може призвести до великого значення функції (U1Д)  і навіть переповнення розрядної сітки ЕОМ (мал. 8). На мал.9 показані для довільної функції (Х) випадки зациклювання ітерацій і розбіжності методу Ньютона. Проте з геометричних інтерпретацій видно, що якщо початкове наближення вибране близько до точного рішення, то метод сходиться завжди.

Мал.8

Мал. 9

Для усунення невиправданого зростання (U1Д) і переповнення розрядної сітки ЕОМ у разі експонентних нелінійностей існує декілька способів [l, 2]:

а) введення обмежень на зміну напруги і струму діодів:

макс , макс 

б) лінеаризація діодних характеристик після Uд макс,  тобтозастосування співвідношень:

           Iо•(eUд/т -1)                            при Uд макс = Uдм             

 Iд=

           Iо•(eUдм/т -1)•(1+(Uд-Uдм)/ Т)        при  Uд > Uдм

в) використання допоміжних співвідношень - визначення поправки, наприклад, при  U >0 за формулою:

UкД =  Т •Ln(1+UкД / Т)

де UкД - поправка, обчислена за звичайною ітераційною схемою Ньютона. Ці ідеї можуть бути перенесені на інші класи нелінійних функцій.

МЕТОДИ РОЗВЯЗАННЯ СИСТЕМ ЛІНІЙНИХ АЛГЕБРАЇЧНИХ РІВНЯНЬ, ЩО ЗАСТОСОВУЮТЬСЯ ДЛЯ АНАЛІЗУ СХЕМ

 

Під час розрахунку статичного режиму методом Ньютона (див. (33)) виникає необхідність розв’язання системи лінійних алгебраїчних рівнянь на кожній ітерації. Чисельні методи розвязання систем лінійних рівнянь поділяються на дві групи:

1) точні (прямі) методи;

2) ітераційні.

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

Ітераційні методи застосовуються, як правило, для ітераційного поліпшення розвязків, що отримуються прямими методами. Ітераційні методи є наближеними. Вони дають розвязок системи як межу послідовних наближень, що обчислюються деяким одноманітним процесом (наприклад, розглянутий метод Ньютона, метод простої ітерації, метод Некрасова, метод Зейделя та інш. [4,5]).

Найпростіший серед точних методів - метод Гауса [4,5]. Він заснований на ідеї виключення невідомих, внаслідок якого задана система рівнянь

W•Х=       або      А•Х=В

тобто

а11•х112•х213•х3+……а1n•хn=b1

а21•х122•х223•х3+……а2n•хn=b2                    (37)

  ………………………………………

аn1•х1+ аn2•х2+ аn3•х3+……аnn•хn=bn

перетворюється у еквівалентну їй систему з верхньою трикутною матрицею, розвязання якої вже не викликає складностей. Метод Гауса може бути реалізований таким чином. Нехай а110, розділимо перше рівняння системи (37) на коефіцієнт а11, який має назву ведучого для першого кроку. Потім помножимо послідовно отримане рівняння на аi1 та (i=2,3,..., n)  і віднімемо його з відповідних рівнянь (i=2,3,..., n)  системи (37). У результаті невідоме x1  буде виключено з усіх рівнянь заданої системи окрім першого, і ми отримаємо систему, еквівалентну (37) вигляду

х112(1)•х213(1)•х3+……а1n(1)•хn=b1(1)

0+а22(1)•х223(1)•х3+……а2n(1)•хn=b2(1)                            (38)

  ………………………………………………

0+аn2(1)•х2n3(1)•х3+……аnn(1)•хn=bn(1)

З цією системою обходимось аналогічно, але без урахування першого рівняння. Таким чином, на другому кроці фактично перетворюваною є система     ( n-1 )-го порядку з матрицею

| а22(1)  а23(1)   …… а2n(1)   |

| а32(1)  а33(1)   …… а3n(1)   |

| ………………………..  |

| аn2(1)   аn3(1)   …… аnn(1) |

і правою частиною b2(1)   b3(1)   …… bn(1) |t

Після другого кроку отримуємо систему, в якій х2 буде виключено з усіх рівнянь, окрім першого і другого. Продовжуючи описаний процес, після n -го кроку прийдемо до системи, еквівалентної (37), але з трикутною матрицею

х112(1)•х213(1)•х3+……а1n(1)•хn=b1(1)

0   +     х2  +………………а2n(2)•хn=b2(2)                                (39)

  ………………………………………

…………………………………хn=bn(n)

Перетворення системи (37) у систему (39) називається прямим ходом, а розвязання трикутної системи (39) – зворотнім ходом. Обчислювальні формули цього варіанту методу Гауса, який має назву алгоритма єдиного поділу, мають наступний вигляд:

Прямий хід,   s -й крок (s = 1, 2,……n)

аiк(s)= аiк(s-1)/аss(s-1), bi(s)= bi(s-1)/аss(s-1), i=s,    к=s,s+1,……,n                                 (40)       

аiк(s)=аiк(s-1)-[аsк(s-1)/аss(s-1)]•аis(s-1),  

bi(s)=bi(s-1)-[аis(s-1)/аss(s-1)]•bi(s-1),                  i=s+1, s+1, ………., n    к=s, s+1, ……, n    

                                    

Зворотній хід здійснюється за формулою:

             n               

xi=bi(i)-аiк(i)•xк,           i=n, n-1, ………,1                                            (41)

             к=i+1            

Схема єдиного поділу проста і економна за кількістю арифметичних операцій (вимагає множень - (n3+3•n2+n)/3, додавань - (2•n3+3•n2+5•n)/6, ділень - n ), однак для її застосувнння необхідно, щоб усі ведучі елементи    аss(s-1) (s=1,2,….,n) були відмінні від нуля. Близькість ведучих елементів до нуля може призвести до великої втрати точності обчисленого розвязку. У зв'язку з цим вводяться різні варіанти методу Гауса, наприклад, алгоритм з вибором головних елементів по усій матриці. Порядок виключення невідомих у заданій системі відбувається таким чином. На кожному кроці s (s=1,2,….,n-1) з коефіцієнтів перетворюваної матриці вибирається найбільший за модулем, він  має назву головного елементу s -го кроку. Невідоме при ньому виключається за описаним вище правилом. Дня зручності обчислень перед виключенням цього невідомого роблять перестановку рівнянь і невідомих так, щоб головний елемент зайняв лівий верхній кут перетворюваної матриці. Якщо на s-му кроці найбільший елемент вибирається серед коефіцієнтів s-го стовпчика (рядка), то такий алгоритм називається алгоритмом з вибором головного елемента по стовпчику (рядку). Потрібно відмітити, що процедура обернення матриць шляхом застосуання виключень Гауса вимагає приблизно n3 множень у порівнянні з n3/3 при розв’язанні системи лінійних рівнянь. Тому не пропонується розвязувати рівняння WХ= шляхом обернення W.  

      При розв’язанні систем лінійних рівнянь, у тому числі під час аналізу лінійних схем, широке застосування крім методу Гауса отримали також LU - розкладання, метод Краута, методи відображень і обертання [5, Д2]. Під час розрахунку великих електронних схем матриця W має значну кількість нульових елементів (до 80-90%) тобто сильно розряджена. Врахування цього факту у спеціальних модифікаціях вищезгаданих методів [5, Д2]  дозволяє різко збільшити ефективність розв’язку (зменшити витрати пам'яті і збільшити швидкість розв’язку).

Під час розрахунку статичного режиму ТРП (див. (9)) за умов статики (33) та із застосуванням методу Ньютона знадобиться розвязувати систему лінійних алгебраїчних рівнянь вигляду

-1/R1

к1

=

1(к1)

1/RБ+

+1/R2

1/R3  -

-1/RБ

к2

2(к2,к3)

1/RБ

W33

N'э-

-

-I'к- -

к3

3(к2,к3,

к4,к5)

N'э++

W44

к4

4(к3,к4)

'к+

+N

-N

-'к+

+1/R4

к5

5(к3,к4,

к5)

1/R6

к6

6(к6)

1/R2

1/R7

-1

iкЕ

7(к2,к5,

iкЕ)

де    'к=dккi/dкi ,       'э=dэкi/dкi    

i - індекс відповідного потенціала,

W33= -1/RБ+'э+N'э+'к+I

W44= 1/R5+'э - N'э

Якщо немає необхідності обчислення iЕ, то можна розвязувати систему 6-го порядку після попереднього видалення останнього рядоку і сьомого стовпецю в матриці W. Потенциали на кожній ітерації визначаються з

к+!i =кi - кi,   к=0,1,2,……..

Якщо [ik ] <  и [(i k)] < ,  то  ik+1 = *i      тобто шуканаму розв’язку.

Для моделей ММС-Ст1 і ММС-Ст2 може бути побудована ітераційна схема (33) подібно до того, як це було зроблено для моделі ТРП. Наприклад, матриця Якобі від (27) має вигляд

W(iЕ ,)=|АЕ|+|АR|•|G|•|АtR|+|АH|•'(|АtH|•||)•|АtH|+

+|АД|•'Н(|АtД|•||)•|АtД|          (42)

Щоб розрахувати методом Ньютона статичний режим заданої схеми з класу електронних схем, представленого множиною гілок (2)-(3), необхідно за допомогою заданої топологічної інформації (матрацям |А| або |П|) сформувати не тільки моделі (27-29) або (31), але й матрицю Якобі. У ЕОМ це здійснюється автоматично за допомогою спеціальних алгоритмів формування моделі і Якобіана. Якщо під час розвязання лінійної системи застосовується метод розріджених матриць, то Якобіан подається не у формі матриці [n x n], а у вигляді звязаних векторів (множин) [5].

Окрім методу Ньютона під час розрахунку статичного режиму, тобто при розвязанні ММС-Ст1 і ММС-Ст2, можна застосувати інші методи:

- метод Ньютона-Канторовича, що вимагає обчислення Якобіана тільки один раз;

- метод найшвидшого спуску, який менш критичний, ніж метод Ньютона, до вибору початкового наближення, але з меншою швидкістю збіжності;

- метод Матвеєва – комбінацію зі згаданих вище;

- методи високої швидкості збіжності [2,4,5].

На основі ММС-Ст розраховується ряд якісних показників, наприклад, споживана потужність на постійному струмі- Ро, як сума потужностей, споживаних у кожній гілці, передавальна характеристика Uвых=(Uвх)  при лінійному законі зміни Uвх та інші. По закінченню розрахунку статичного режиму знаходять також параметри малосигнальних моделей приладів  Наприклад, для транзистора визначається

                                        gK= I'К•(U*КБ)  IК/U*КБ |

                    | UКБ=U*КБ

rK=1/gK,    CK=CK(UКБ)

gЭ= I'Э•(U*БЭ)  IЭ/U*БЭ    |

                    | UБЭ=U*БЭ

rЭ=1/gЭ,   CЭ=CЭ(UБЭ)

U*КБ, U*БЭ - значення напруг у статичному режимі.

МАТЕМАТИЧНА МОДЕЛЬ ЕЛЕКТРОННОЇ СХЕМИ У ДИНАМІЧНОМУ РЕЖИМІ ПРИ МАЛОМУ СИГНАЛІ.

ЧАСОВА ОБЛАСТЬ

При малому вхідному сигналі нелінійні функції, що описують властивості приладів, в околі статичного режиму можна вважати лінійними і представляти параметрами, розрахованими (або виміряними) у статичному режимі. Тоді на еквівалентній схемі транзистора (мал. 2) діоди будуть замінені резисторами rK и rЭ ,  нелінійні ємності - постійними значеннями. Математична модель ТРП (9) набуває вигляду лінійних алгебро-диференціальних рівнянь, а її узагальнений матричний аналог (14) -ММС-ДР1

Е|•|iЕ|+|АR|•|G|•|АtR|•||+|АH|•gН•|АtH|•||+|АД|•gД•|АtД|•||+

+|АСЛ|•|СЛ|•d(|АtCЛ|•||)/dt+|АСН|•|СНd(|АtCН|,||)/dt = 0  (43) 

gН=|1/rк  1/rэ|,    СН=|Ск  Сэ|

або після перестановок:

d(Асл, )/dt)

|Асл Cл|Асн Cн |АЕ|

d(Асн, )/dt)

=

[ARGRAtR+AHgHAtH+ AДgДAtД] +

+ AR1G1•UВХ                                                 (44)

iЕ

и, нарешті, у скороченому вигляді без рівняння для iЕ

A•d(AtC, )/dt=AG+ AUu(t)        (45)

Відповідно ММС-ДР2 (див.(23)) для малого сигналу після алгебраїчних перетворень з урахуванням розбиття матриць tRХ і tIХ на субматриці запишеться як:

С•duc/dt=|C|-1•[(CRXGХtCRX+CIHgHtCIX+CIДgHtCIX)•uc+ +(CRXGХtRBRX+CIHgHtRBIX+CIДgHtRBIX)•uRB+

+(CRXGХtERX+CIHgHtEIX+CIДgHtEIX)•EB]                (46)

або більш скорочено

 dX/dt=A1X +A2X+B1u(t),    X=uC, Xл=uRB, u(t)=EB             (47)

де A1, A2, B1  - матриці у квадратних дужках, помножені на |C|-1.

Використовуючи подібні перетворення з (23) і (24), маємо

XЛ1X2u(t)           (48)

Підставляючи (48) у (47) остаточно отримаємо

dX/dt=(A1+A2Д2)•X+(B1+A2Д2)•u(t),    

або

dX/dt=A•X+B•u(t)=( X, t )                      (49)

Це рівняння подане явно відносно похідних dX/dt,  а рівняння (45) представлене у неявній формі відносно похідних, порядок і вигляд систем рівнянь також різний. Алгоритм формування рівнянь (45) більш простий, бо не вимагає вибору дерева на графі і такої кількості операцій з матрицями, як при отриманні (48).

Задача розрахунку динамічного режиму математично формулюється як задача Коши, яка полягає у тому, що шукається розв’язок X(t) (або (t))  рівняння (49) (або (45)), який би задовольняв задану початкову умову

X(tО)=XО=X*((tО)=*) 

на інтервалі tО  t  tК,   tК - кінець інтервала.

Умови існування і єдиності розвязку поставленої задачі Коши будемо вважати виконаними.

Під час розв’язку чисельними (кінцево-різностними) методами шукана функція шукається в окремих точках інтервалу  [tО, tК],  t0,t1,…,tn, …,tN, які мають назву вузлів, у вигляді таблиці значень X0,X1,…,Xn,…,XN, що наближено дорівнюють значенням X(t1),X(t2),…,X(tn),…,X(tN) точного розв’язку X(t).  Відстань між вузлами h=t=tn+1-tn має назву крока    інтегрування і може бути задана або перед початком обчислень (інтегрування з постійним кроком), або визначатися у процесі обчислень (інтегрування з автоматичним вибором кроку).

Більшість чисельних методів розвязання задачі Коши, що розглядається, можна подати у вигляді [l, 2,4.5, Д4]

Xn+1=F(Xn-q, Xn-q+1,…,Xn,Xn+1,…, Xn+s)           (50)

де F - деяка відома функція вказаних аргументів, що визначається способом побудови методу і залежить від інтегрованого рівняння і вибраної сітки

tО < t1 <….. <tN=T

При q=0,    0 s 1   такі обчислювальні алгоритми називають однокроковими, а q 1 або s 1 – багатокроковими. Як однокрокові, так і багатокрокові методи вигляду (50) називають явними у випадку s=0  і неявними при s=1. У випадку, коли s=1, багатокрокові алгоритми називають методами із забіганням вперед. Багатокрокові методи вимагають застосуання спеціальних обчислювальних алгоритмів для знаходження перших q- значень X1,X2,..,Xq

наближеного розв’язку і останніх його s-1   значень XN-S+2,X N-S+3,…,XN.

У зв'язку з наближеним характером розвязку на кожному кроці чисельного інтегрування виникає так звана локальна помилка:

к =|| XК, X(tК)||

між точним і отриманим розвязком. Локальна помилка складається завжди з двох компонент:

1) помилки методу (-або відкидання) -мк 

2) помилки округлення

Перша залежить від виду чисельного алгоритму, що використовується при    обчисленні XК , а друга зумовлена кінцевою довжиною машинного слова при реалізації алгоритму на ЕОМ. Як помилка методу, так і помилка округлення нагромаджуються із збільшенням числа кроків. Деякі методи (чисельно нестійкі) сприяють посиленню локальної помилки методу і помилки округлення на кожному кроці так, що через деякий час помилка, що зросла, може переважати над самим розвязком. Стійкість же методу гарантує, що локальні помилки не посилюються, а залишаються обмеженими для досить малої величини кроку h при h. Методи мають різну область стійкості і, відповідно, різну величину максимально-можливого кроку інтегрування. На результат розвязку впливає також точність завдання початкових умов. На практиці для розрахунків доцільно застосовувати методи, що сходяться.  Метод вважається таким, що сходиться в тому випадку, коли для розв’язання задачі Коши, яка має єдиний розв’язок X(t), обчислений розв’язок X(tК) при дії вказаних чинників однозначно збігається до X(t) на інтервалі tО  t  T при  t і h=T/n.

Найпростішими з методів чисельного інтегрування є явний і неявний методи Ейлера. Вони є першими у ряду як однокрокових, так і багатокрокових методів. Згідно з явним методом Ейлера розвязок визначається за формулою:

 

Xn+1= Xn + h(Xn, tn ) + rn+1          (51)

де    Xn+1=X(tn+h),    tn= to+n•h,

rn+1- помилка методу (або помилка відкидання),

rn+1=(h2/2)•''(tn+•h, Xn)  0 <  < 1           (52)

або виражена через кінцеві різниці значення Xi

rn+1=(h2/2)•(2•Xn)/h2= (Xn+1-2•Xn+Xn-1)/2          (53)

Показник при h (у нашому випадку два) характеризує порядок точності методу.

Геометрично метод Ейлера означає заміну інтегральної кривої, що являє собою точний розв’язок X(t) кусково-ламаною лінією, ділянки якої паралельні дотичній до X(t) у вузлах tn. Тому метод Ейлера називають інакше методом ламаних або дотичних.

Якщо помилка задана, тобто rn+1= , то крок інтегрування  може бути вибраний виходячи з похибки методу

h [2•/(''(tn,Xn))]0.5                              (54)

Значно більші обмеження на крок інтегрування явного методу Ейлера накладаються з умов стійкої поведінки у процесі обчислень [2,6]. Так при інтегруванні стійкої системи лінійних диференціальних рівнянь (50), що має, наприклад, р різних дійсних від’ємних власних чисел матриці А=|1,2,…,p|t, розв’язок прямує до нуля, тобто LimnXn=0 тільки за умови виконання 

|1+h•i|t <1,        i=1,2,……,p

Враховуючи  i  отримуємо верхню межу для величини кроку інтегрування

h < 2/МАКС           (55)

де   МАКС = МAX[|1|,|2|,|3|,….,|p1|,]

Відомо, що постійні часу лінійної схеми, що описується системою (50), пов'язані з власними значеннями |А| наступним чином:

=-1/i.

Тому h < 2•МИН, де МИН =МIN[|1|,|2|,|3|,….,|p1|,]. Якщо у системі рівнянь є великий розкид власних значень i (кажуть, жорстка система) або у схемі великий розкид постійних часу, то потрібно інтегрувати з малим кроком навіть на ділянці, де розвязок змінюється повільно (велике ). Будь-яка спроба збільшення кроку негайно призводить до різкого зростання похибки ( "вибуху" похибки).

Неявний метод Ейлера       Xn+1= Xn + h(Xn+1, tn+1 )                       (56) 

має похиьку метода                       rn+1= -0.5•h2''(Xn, tn),                            (57)

але на відміну від явного методу Ейлера його властивості стійкості не накладають яких-небудь обмежень на крок.   Дійсно, з умови стійкості

1/|1-h•i|   <  1

і оскільки  i<0, то наближений розвязок стійкий для всіх h>0. Крок інтегрування може бути вибраний виходячи тільки з міркувань точності, тобто по формулі (57). Застосування неявного методу Ейлера для розв’язку жорстких систем рівнянь (з великим розкидом постійних часу) дозволяє отримати виграш у кількості кроків h  і тим більший, чим більше розкид постійних часу. Однак використання (56) пов'язане з розвязанням на кожному кроці інтегрування системи нелінійних рівнянь для знаходження вектора Xn+1, якщо  (Xn) - нелінійна.

Застосування (56) для інтегрування (49) призводить до розвязку на кожному кроці системи

|1/h - A |•Xn+1=(1/h)•Xn+1 + B•u(t n+1)                (58)

а для розвязку  (45) - системи лінійних рівнянь

|A/h - AG |• n+1= (A/h)• n + Auu(t n+1)         (59)

При явному методі Ейлера отримуємо, відповідно,системи

 Xn+1=|1+h•A |•Xn + h•B•u(t n1)          (60)

(A/h)• n+1=|A/h + AG |• n + Auu(t n)          (61)

Розвязання лінійних систем може бути зроблене методом Гауса або іншим методом розвязку систем алгебраїчних рівнянь. Треба зазначити, що при інтегруванні жорстких систем, в лінійних рівняннях можуть виникати погано зумовлені матриці, що висуває додаткові вимоги до методів розв’язку лінійних систем і обмежує крок h.

Методи Ейлера, як правило, при аналізі динамічного режиму складних електронних схем не завжди забезпечують необхідну точність і ефективність розвязання. Тому застосовують методи більш складні, більшого порядку точності і у разі розв’язання жорстких систем, а такими у більшості випадків є математичні моделі електронних схем, - жорстко-стійкі:

  •  метод Шихмана,
  •  метод Гира,
  •  ФДН,
  •  неявні методи Рунге-Кутта [1,2,5].

У випадку коливального характеру розвязку перевага надається методу трапецій і неявним методам Рунге-Кутта.

Маючи розв’язки  (58)-(61) і рівняння (13), (23)-(24) або (49) можна визначити майже усі струми і напруги у схемі, що аналізується, і необхідні якісні показники: коефіцієнти передачі, вхідні і вихідні опори та інше.

МЕТОДИ І АЛГОРИТМИ РОЗРАХУНКУ ДИНАМІЧНОГО   РЕЖИМУ ПРИ ВЕЛИКОМУ    СИГНАЛІ.

ЧАСОВА ОБЛАСТЬ.

Математичні моделі для динамічного нелінійного режиму, детально розглянуті у розділі "Математична модель електронної схеми у динамічному режимі при великому сигналі. Часова область", (9), (13), (14), (23)-(25), можна подати як  (10) і (26) у двох характерних узагальнених формах

Fд (X, X, t)=0,   X(to)=X*           (62)

 і

X= д (x, t),               (63)

де , д - нелінійні оператори,

X* - результат розрахунку статичного режима,  X|=|XД,XН,XЛ|t,

тобто у вигляді нелінійних алгебро-диференціальних рівнянь і у нормальній формі звичайних диференціальних рівнянь. Розвязок цих рівнянь може бути знайдений при використанні явних і неявних методів чисельного інтегрування. Застосовуючи явний метод Ейлера до (63) знайдемо розвязок у вузлах t0, t1, t2,….., tn, інтервала [t0,….. tК]  по співвідношенню

 Xn+1= Xn+ h•д(Xn, tn )          n=0,1,……N               (64)

З рівняння (56) для похідної Х отримаємо

Xn+1=(Xn+1, tn+1 )=(Xn+1 - Xn)/h             (65)

Введення з (65) у (62) і (63) дозволяє останні звести до кінцевих нелінійних алгебраїчних рівнянь у вузлах різносної сітки

Fд (Xn+1, (Xn+1 - Xn)/h, tn+1)=0          (66)

(Xn+1 - Xn)/h=д(Xn+1, tn+1 )                     (67)

Для розвязку (66) і (67) застосуємо метод Ньютона. У результаті отримаємо алгоритми розрахунку динамічного режиму у вигляді

|Fд(Xкn+1)/Xкn+1|•Xкn+1= - Fд(Xкn+1, (Xn+1 - Xn)/h, tn+1)      (68)

Xк+1n+1=Xкn+1 + Xкn+1,     к=0,1,2,….,

і

 |д(Xкn+1)/Xкn+1-1/h|•Xкn+1=д(Xкn+1, tn+1) - Xкn+1/h + Xn,         (69)

Xк+1n+1=Xкn+1 + Xкn+1

Таким чином, розрахунок динамічного режиму зводиться до розвязання к-раз на кожному кроці чисельного інтегрування систем лінійних алгебраїчних рівнянь, наприклад, розглянутим раніше методом Гауса. Умови збіжності методу Ньютона вимагають досить доброго початкового наближення, точність якого і визначає число ітерацій у (68) і (69). Для його отримання рекомендують [2,5]один з явних методів (наприклад, метод Ейлера). При цьому керуються тим, що характеристики стійкості явного методу не мають істотного значення при його одноразовому застосуванні. Треба зазначити, що з умов збіжності методу Ньютона витікають додаткові обмеження на межі зміни кроку чисельного інтегрування h. При побудові чисельних алгоритмів треба також враховувати особливості інтегруючих систем рівнянь (наприклад, жорсткість), які можуть призвести до поганої обумовленості матриці Якобі і труднощів під час розв’язання систем лінійних рівнянь. При сильно розрідженій матриці Якобі застосовують методи розвязання систем з розрідженими матрицями  [4,5].

Компоненти вектора Fд   і матриці Якобі Fд/X для ММС  ТРП вигляду (9) при підстановці (69) будуть:

  

1=(С1/h-1/R1)•1(n+1)-2/h)•2(n+1)-uвх(tn+1)/R1+(С1/h)•1n+(С2/h)•2n

2=(С1/h)•1(n+1)-1/h+1/RБ+1/R2)•2(n+1)+(1/R3-1/RБ+1/R2)•3(n+1)-

-E2/R2+(С1/h)•1n+(С2/h)•2n

3=…………………………………………………………………………         (70)

4= -[ СЭ(-3(n+1)+4(n+1)) / h ]•3(n+1)+[ СЭ(-3(n+1)+4(n+1)) / h ]•4(n+1)-

-Iк(-3(n+1)+4(n+1)) +э(-3(n+1)+4(n+1)+5(n+1)/R5 --

-СЭ(-3(n+1)+4(n+1)) / h ]•3(n+1)+СЭ(-3(n+1)+4(n+1)) / h ]•4(n+1)

Fд 5=………………………………………………………………………………..

Fд 7=2(n+1)/R2 -- 5(n+1)/R4 - iE(n+1)-- (1/R2 +1/R4)•E                                                  (71)

1

2

3

4

5

6

7

1

С1/h-1/R1

1/h

……………………………..

..

……

.

….

2

С1/h

1/h+1/RБ+1/R2

…………………………….

……

.

….

Fд(Xкn+1)/Xкn+1=

3

………..

…………….…

…………………………….

……

.

….

4

………

………………

(С'Э/h)•3(n+1)+(С'Э/h)+I'к--'э

..

.

….

5

………..

……………….

……………………………………

….

.

….

6

……….

………………………………….

..

.

….

7

………….

1/R2

…………………………………

-1/R4

-1

 Тут 'к, 'э, С'э - похідні від нелінійних функцій по відповідним змінним.

Узагальнені моделі (14) і (25) таким же шляхом, як було зроблено для (9), можна подати і вигляді (68) і, відповідно, (69).

На основі системи (14) вектор Fд  і матриця Якобі Fд/ запишуться у формі:

Fд= AСН•[CН(AtСН,n+1)/h]•AtСНn+1 +

+ AСЛ•(CЛ/h)•AtСЛn+1 +

+ AНн(AtН, n+1) +

+ AНн(AtД,n+1)+

+ AR•GR•AtRn+1 +            (72)

+AEiК +

+AСН•[CН(A*СН, n)/h]•AtСНn +

+ AСЛ•(CЛ/h)•AtСЛn

 Fд/n+1 = AСН•[C 'Н(AtСН, n+1)•AtСН /h]•AtСНn+1 +

+ AСН•[CН(AtСН, n+1) /h]•AtСН +  

+ AСЛ•(CЛ/h)•AtСЛ +

+ AН 'н(AtН, n+1)•AНt +             (73)

+ AН 'н(AtД, n+1)•AtД+

+ AR•GR•AtR+

+ AE    

У процесі реалізації на ЕОМ задачі розрахунку динамічного режиму для електронної схеми довільної структури вектор Fд  і Якобіан Fд/n+1     формуються по заданій топології схеми автоматично за допомогою спеціальних алгоритмів формування.

Алгоритм аналізу динамічного режиму, що застосовує ММС вигляду (63) або (25), і явні методи чисельного інтегрування не вимагають розв’язання систем нелінійних алгебраїчних рівнянь і, здавалося б, більш ефективні у витратах часу і пам'яті на ЕОМ, ніж (68) або (69). Але, насправді, через обмеження на крок, з умов стійкості явних методів, цей алгоритм, особливо при інтегруванні жорстких систем диференціальних рівнянь (з великим розкидом постійних часу), а такими майже завжди є рівняння електронних схем, виявляється менш ефективним.

МАТЕМАТИЧНА МОДЕЛЬ ЕЛЕКТРОННОЇ СХЕМИ У ДИНАМІЧНОМУ РЕЖИМІ ПРИ МАЛОМУ СИГНАЛІ.

ЧАСТОТНА ОБЛАСТЬ.

Математичні моделі для частотної області вотримують на основі моделей (43), (45), (49) для часової області за допомогою перетворення Лапласа (оператор d/dt замінюється на p або j , а змінні (t)(p), X(t)X(p), u(t)u(p) і урахуванні того, що опір джерела живлення на змінному струмі дорівнює нулю, тобто Е=0).

Тоді рівняння (45), наприклад, набуває вигляду

|A•p•AtC|•(p)=|AG|•(p)+|AU|•u(p)         (74)

або

|A•p•AtC  - AG|•(p)=|AU|•u(p)            (75)

а з (47)-(49) можна отримати

|p - A|•Х(p)=В•u(р)              (76)

Хл(р)=[Д1•|p-A|-1•B+Д2]•u(p)             (77)

Якісні показники або функції схеми (наприклад, коефіцієнти передачі) визначаються з (75)-(77) таким чином

  K(p)=[(p)/u(p)]•[ Au/|A•p•AtC  - AG |]                            (78)

  K1(p)=[X(p)/u(p)]•[В/(р-А)]                  (79)

  K2(p)=[Xл(p)/u(p)]•[ Д1•В/(р-А)+Д2]                  (80)

Математична модель може бути складена і безпосередньо по еквівалентній схемі змінного струму. Наприклад, для ТРП, який розгядається, еквівалентна схема, граф схеми і матриця інциденцій А подані нижче.

Мал.10

Еквівалентна схема ТРП для режиму малого сигналу в частотній області

Мал. 11

Граф схеми ТРП для режиму малого сигналу в частотній області

Матриця інциденцій |А|

1

2

3

4

5

6

7

8

9

10

1

-1

1

2

-1

-1

1

A

=

3

-1

-1

1

4

1

-1

5

1

-1

1

6

-1

-1

Кожна k-та гілка схеми подається як узагальнена

iК =bК•(uкк) -Jк

де bК - провідність гілки,   в загальному  bК=1/R+p•C+1/(p•L)  

Ек, Jк – джерела струму і напруги в гілці, iК – струм у гілці

uк- напруга на bК .

Вузлове рівняння вигляду (63) для даной схеми буде

|Y|•||=|I|        (81)

де    |Y|=|A|•|YB|•|A|t,   |I|=|A| • (|J|-|YB|•|E|)

|YB|- матриця провідностей гілок,

|E| и |J| - вектори джерел напруги і струму гілок,

||=|1, 2, …., 6,|t - вектор потенціалів вузлів.

Матриця провідностей вузлів |Y|

1/R1+pC1

-pC1

-pC1

pC1+1/R2+

+1/R3+1/RБ

-1/RБ

-1/RБ

p•(CЭ+CК)+

+1/rЭ+1/RБ+

+1/rК

-1/rЭ+pC3

Y

=

-1/rЭ+pC3

p•(CЭ+C2)+

+1/rЭ+1/R5

pC4+1/rK+

+1/R4

-pC3

-pC3

p•(C3+C4)+

+1/R6

Вектор еквівалентних вузлових джерел струму

|I|=|-uвх/R1, 0, Iэ, 0, Iэ, 0|t.

Напруга вузлів і гілок                    

|u|=|A|t•||,  ||=|1, 2, …., 6,|t,   |u|=|u1, u2, …., u10,|t .

Коефіцієнт передачі за напругою з вихода схеми на вхід визначиться як

Кu=6/uвх             (82) 

а коефіцієнт передачі за потужністю – з виразу

 Кu=[26(Z11+R1)] / [ u2вхZ66]          (83)

де Z11 і  Z66  - відповідні елементи Y-1

Пошук розвязку систем рівнянь (75), (76), (81) і визначення якісних показників (див. (78)-(83)) пов'язані з розв’язком систем лінійних рівнянь (або зоберненням матриць) з комплексними коефіцієнтами. Описані раніше алгоритми для розв’язку систем лінійних рівнянь можуть бути застосовані і у цьому випадку з урахуванням того, що усі величини комплексні. Повна кількість операцій множення тут у 4 рази більше, ніж з дійсними коефіцієнтами. Це видно з виразу при множенні двох комплексних чисел Z1,  Z2

Z1Z2=(x1+jy1)•(x2+jy2)= (x1•x2 - y1•y2) + j(x1•y2 + y1•x2)

МАТЕМАТИЧНА МОДЕЛЬ ЕЛЕКТРОННОЇ СХЕМИ У ДИНАМІЧНОМУ РЕЖИМІ ПРИ ВЕЛИКОМУ СИГНАЛІ.

ЧАСТОТНА ОБЛАСТЬ.

При впливі на схему великого гармонічного сигналу  для складання математичної моделі схеми в частотній області застосовують однократне перетворення Фурє.

 Якщо на нелінійний елемент з характеристикою

(uн1,uн2,…,uнz)

впливає сигнал великої потужності, то розкладання в ряд Фурє має вигляд:

                        N               

        н(t) =(1/2) •  Фn ej(n•)•t       (84)

                     n=-N             

де

                  2

                   

Фn=(1/2)•( uн1,uн2,…, uнz) • е-j(n1)•d1 ,        1=1•t,   (85)

                   

           0        

                        N               

uнq(t)=(1/4)  Unнq ej(n•1)•t  q=1,2,……,z 1=1•t,

                      n=-N

Фn и Unнq - комплексні амплітуди гармонік;

n- номер гармоніки 1=n•;

N- кількість гармонік, що використовуються;

Для спрощення запису введемо  позначення   wn= еj(n1)    (86)

              N                                           N           w-n= е-j(n1)

а  замість       будем использовать   

           n=-N                 n

врахуємо також 

N

  u(t)=(1/4)•  Unнq ej(n•1)•t     (87)

n

      N

  iвх(t)=(1/4)•  Inвх • w-n     (88)

 n         N

при переході в частотну область   d/dt  j(r•1)

             r

Для нелінійного конденсатора   iсн(t)=С(uс)•duс/dt   (89)

Комплексні амплітуди визначаються з співвідношення

                     2

                                         n           

Icn=(j/42)•С(uсн)•wn[r Urc•w-r]•d1      (90)

                                         r

              0         

напруга на нелінійному конденсаторі

                     N

uсн(t)=(1/4)  Ucn  w-n        (91)

                  n

струм через нелінійний конденсатор

                     N

iсн(t)=(1/4)  Icn  w-n         (92)

                  n

струм через нелінійний резистор

                    N

iн(t)=(1/4)  Inн  w-n         (93)

                 n

де

                   2

                              N           

Inн=(1/2)• iн(t)( Unн•w-n)• wn d1       (94)

                              n

            0         

Для лінійних компонент запишемо

лінійна ємність

                     N

iсл(t)=(1/4)  Inсл  w-n        (95)

                  n

де              N

     Inсл =j(r1) •C•Urc w-r       (96)

                 r

лінійна провідність

                   N

 iG(t)=(1/4)GUnG w-n        (97)

                n

З урахуванням всіх вище зазначених висновків в матричній формі отримаємо систему в області комплексних амплітуд, використовуючи одномірне перетворення Фур’є

                N                           N                       N                        N

(|Асн|/4)•  Inсн +(|Ас|/4)• Inс+(|Ас|/4)• Inн+(|АR|/4)•G|Ас|tUn=

           n                            n                        n                        n

                                         N

                       = (1/4) •Inвх       (98)

                                         n

де

                        2

                                               N                       N

Incн=(j/42)•С• (|Асн|/4)•[Un•w-n]•wn[r|Асн|tUr•wr]•d1  (99)

                                               n                              &nbs