Расчет переходного процесса в RLC-цепи методом Рунге-Кутта, страница 9

Текст программы для расчета данных на языке QBASIC

DECLARE SUB tables (c!,r!)

CLS

FOR c = .001 TO .0022 STEP .0004

r=7

CALL tables(c,r)

NEXT c

FOR c = .001 TO .0022 STEP .0004

r=70

CALL tables(c,r)

NEXT c

END

SUB tables (c,r)

f$ = "TABLES.txt"

OPEN f$ FOR APPEND AS #100

REM Исходные данные

 E = 100

 L = 1.5

 time = 2 * 2 * 3.14 * SQR(L * c)

 h = .1 * 2 * 3.14 * SQR(.001 * L)

 UC1 = 0

 IC1 = 0

 UR = 0

 UC = 0

 IC = 0

 UL = E

PRINT #100,

PRINT #100, USING "R =#.#Ом  C = #.####Ф"; r; c

PRINT #100,

PRINT #100, "  t      Uc      Ic       Ul       Ur      K1u       K1i     K2u       K2i     K3u       K3i     di     du   "

FOR j = 0 TO time STEP h

UC = UC1 + (k1u + 4 * k2u + k3u) / 6

IC = IC1 + (k1i + 4 * k2i + k3i) / 6

UL = L * (IC - IC1) / h

IF j = 0 THEN UL = E

UR = r * IC

k1u = h * IC / c

k1i = h * (E - r * IC - UC) / L

k2u = h * (IC + k1i / 2) / c

k2i = h * (E - r * (IC + k1i / 2) - (UC + k1u / 2)) / L

k3u = h * (IC - k1i + k2i * 2) / c

k3i = h * (E - r * (IC - k1i + k2i * 2) - (UC - k1u + k2u * 2)) / L

di = (k1i + 4 * k2i + k3i) / 6

du = (k1u + 4 * k2u + k3u) / 6

IC1 = IC

UC1 = UC

PRINT #100, USING "#.###  ###.###  ##.###  ####.###  ####.###  ####.###  ####.###  ####.###  ####.###  ####.###  ####.###  ###.###  ###.###"; j; UC; IC; UL; UR; k1u; k1i; k2u; k2i; k3u; k3i; di; du

NEXT j

CLOSE #100

END SUB

Текст программы для построения графиков на языке QBASIC

CLS

E = 100

L = 1.5

DIM C(3), R(2) AS SINGLE

C(1) = .001

C(2) = .0014

C(3) = .0018

R(1) = 7

R(2) = 70

h = .00024

p = 3

w = 2

'строим графики зависимостей параметров цепи от времени

SCREEN 12

FOR w = 1 TO 2

tmax = 2 * 2 * 3.141592 * SQR(.0018 * 1.5)

n = 3

'окно для рисования графиков

VIEW (10, 10)-(630, 470), 0, 1

'разрешение окна

WINDOW (0, 0)-(620, 460)

'строим оси декартовых координат

COLOR 3

LINE (20, 340)-(280, 340), 1

LINE (20, 260)-(20, 440), 1

LINE (20, 40)-(280, 40), 1

LINE (20, 20)-(20, 250), 1

LINE (320, 340)-(580, 340), 1

LINE (320, 300)-(320, 440), 1

LINE (320, 125)-(580, 125), 1

LINE (320, 20)-(320, 260), 1

'подписываем оси

LOCATE 4, 19: PRINT "Ic(t)"

LOCATE 17, 19: PRINT "Uc(t)"

LOCATE 4, 55: PRINT "Ur(t)"

LOCATE 17, 55: PRINT "Ul(t)"

LOCATE 9, 38: PRINT "t"

LOCATE 3, 3: PRINT "Ic"

LOCATE 28, 38: PRINT "t"

LOCATE 16, 3: PRINT "Uc"

LOCATE 9, 75: PRINT "t"

LOCATE 3, 40: PRINT "Ur"

LOCATE 22, 75: PRINT "t"

LOCATE 14, 40: PRINT "Ul"

LOCATE 9, 3: PRINT "0": LOCATE 28, 3: PRINT "0": LOCATE 9, 41: PRINT "0"

LOCATE 22, 41: PRINT "0"

'задаем цикл изменения параметра С

FOR p = 1 TO 3

'задаем значения коэффициентов Рунге-Кутта и параметров цепи при t=0

i = 0: Uc = 0: Ur = 0: Ul = E

t = 0

'запускаем цикл изменения времени

DO

t = t + h

'рассчитываем значения коэффициентов Рунге-Кутта

ip = i

i = i + (k1i + 4 * k2i + k3i) / 6

Uc = Uc + (k1u + 4 * k2u + k3u) / 6

Ur = i * R(w)

Ul = (i - ip) * L / (h)

k1i = (E - i * R(w) - Uc) * h / L

k1u = i * h / C(p)