라플라스 변환을 사용하여 RLC 회로의 미분 방정식 풀기

라이브 스크립트 열기

Symbolic Math Toolbox™에서 다음 워크플로로 라플라스 변환을 사용하여 RLC 회로의 미분 방정식을 풉니다. 라플라스 변환의 간단한 예제는 laplace 및 ilaplace 항목을 참조하십시오.

정의: 라플라스 변환

함수 $f (t)$ 의 라플라스 변환은 다음과 같습니다.

$F (s) = \int_{0}^{\infty} f (t) e^{- ts} dt .$

개념: 기호 워크플로 사용하기

기호 워크플로에서는 계산을 수치 형식이 아닌 자연적인 기호 형식으로 유지합니다. 이 접근 방식은 해의 속성을 이해하고 정확한 기호 값을 사용하게 해줍니다. 수치 결과가 필요할 때나 기호적으로 진행할 수 없을 때만 기호 변수에 숫자를 대입합니다. 자세한 내용은 수치 연산방식 또는 기호 연산방식 선택하기 항목을 참조하십시오. 일반적으로 단계는 다음과 같습니다.

방정식을 선언합니다.
방정식을 풉니다.
값을 대입합니다.
결과를 플로팅합니다.
결과를 분석합니다.

워크플로: 라플라스 변환을 사용하여 RLC 회로 풀기

방정식 선언하기

라플라스 변환을 사용하여 초기 조건이 있는 미분 방정식을 풀 수 있습니다. 예를 들어, 다음 회로와 같은 저항-인덕터-커패시터(RLC) 회로를 풀 수 있습니다.

저항(단위: 옴): $R_{1}, R_{2}, R_{3}$
전류(단위: 암페어): $I_{1}, I_{2}, I_{3}$
인덕턴스(단위: 헨리): $L$
커패시턴스(단위: 패럿): $C$
AC 전압원(단위: 볼트): $E (t)$
커패시터 전하(단위: 쿨롬): $Q (t)$

키르히호프의 전압 및 전류 법칙을 적용하여 다음 방정식을 얻습니다.

$\begin{array}{l} I_{1} = I_{2} + I_{3} \\ L \frac{{dI}_{1}}{dt} + I_{1} R_{1} + I_{2} R_{2} = 0 \\ E (t) + I_{2} R_{2} - \frac{Q}{C} - I_{3} R_{3} = 0 \end{array}$

위의 방정식에 관계식 $I_{3} = d Q / d t$ (커패시터 충전 속도)를 대입하여 RLC 회로의 미분 방정식을 얻습니다.

$\begin{array}{l} \frac{{dI}_{1}}{dt} - \frac{R_{2}}{L} \frac{dQ}{dt} = - \frac{R_{1} + R_{2}}{L} I_{1} \\ \frac{dQ}{dt} = \frac{1}{R_{2} + R_{3}} (E (t) - \frac{Q}{C}) + \frac{R_{2}}{R_{2} + R_{3}} I_{1} \end{array}$

변수를 선언합니다. 물리량이 양수 값을 가지므로, 변수에 이에 해당하는 가정을 설정합니다. $E (t)$ 를 1V의 교류 전압이라고 하겠습니다.

syms L C I1(t) Q(t) s
R = sym('R%d',[1 3]);
assume([t L C R] > 0)
E(t) = 1*sin(t);        % AC voltage = 1 V

미분 방정식을 선언합니다.

dI1 = diff(I1,t);
dQ = diff(Q,t);
eqn1 = dI1 - (R(2)/L)*dQ == -(R(1)+R(2))/L*I1

eqn1(t) = 
 $\frac{\partial}{\partial t} I_{1} (t) - \frac{R_{2} \frac{\partial}{\partial t} Q (t)}{L} = - \frac{I_{1} (t) (R_{1} + R_{2})}{L}$

eqn2 = dQ == (1/(R(2)+R(3))*(E-Q/C)) + R(2)/(R(2)+R(3))*I1

eqn2(t) = 
 $\frac{\partial}{\partial t} Q (t) = \frac{\sin (t) - \frac{Q (t)}{C}}{R_{2} + R_{3}} + \frac{R_{2} I_{1} (t)}{R_{2} + R_{3}}$

방정식 풀기

eqn1과 eqn2의 라플라스 변환을 계산합니다.

eqn1LT = laplace(eqn1,t,s)

eqn1LT = 
 $s laplace (I_{1} (t), t, s) - I_{1} (0) + \frac{R_{2} (Q (0) - s laplace (Q (t), t, s))}{L} = - \frac{(R_{1} + R_{2}) laplace (I_{1} (t), t, s)}{L}$

eqn2LT = laplace(eqn2,t,s)

eqn2LT = 
 $s laplace (Q (t), t, s) - Q (0) = \frac{R_{2} laplace (I_{1} (t), t, s)}{R_{2} + R_{3}} + \frac{\frac{C}{s^{2} + 1} - laplace (Q (t), t, s)}{C (R_{2} + R_{3})}$

함수 solve는 기호 변수에 대해서만 방정식을 풉니다. 따라서 solve를 사용하려면 먼저 laplace(I1(t),t,s)와 laplace(Q(t),t,s)에 변수 I1_LT와 Q_LT를 대입하십시오.

syms I1_LT Q_LT
eqn1LT = subs(eqn1LT,[laplace(I1,t,s) laplace(Q,t,s)],[I1_LT Q_LT])

eqn1LT = 
 $I_{1, LT} s - I_{1} (0) + \frac{R_{2} (Q (0) - Q_{LT} s)}{L} = - \frac{I_{1, LT} (R_{1} + R_{2})}{L}$

eqn2LT = subs(eqn2LT,[laplace(I1,t,s) laplace(Q,t,s)],[I1_LT Q_LT])

eqn2LT = 
 $Q_{LT} s - Q (0) = \frac{I_{1, LT} R_{2}}{R_{2} + R_{3}} - \frac{Q_{LT} - \frac{C}{s^{2} + 1}}{C (R_{2} + R_{3})}$

I1_LT와 Q_LT에 대해 방정식을 풉니다.

eqns = [eqn1LT eqn2LT];
vars = [I1_LT Q_LT];
[I1_LT, Q_LT] = solve(eqns,vars)

I1_LT = 
 $\frac{L I_{1} (0) - R_{2} Q (0) + C R_{2} s + L s^{2} I_{1} (0) - R_{2} s^{2} Q (0) + C L R_{2} s^{3} I_{1} (0) + C L R_{3} s^{3} I_{1} (0) + C L R_{2} s I_{1} (0) + C L R_{3} s I_{1} (0)}{(s^{2} + 1) (R_{1} + R_{2} + L s + C L R_{2} s^{2} + C L R_{3} s^{2} + C R_{1} R_{2} s + C R_{1} R_{3} s + C R_{2} R_{3} s)}$

Q_LT = 
 $\frac{C (R_{1} + R_{2} + L s + L R_{2} I_{1} (0) + R_{1} R_{2} Q (0) + R_{1} R_{3} Q (0) + R_{2} R_{3} Q (0) + L R_{2} s^{2} I_{1} (0) + L R_{2} s^{3} Q (0) + L R_{3} s^{3} Q (0) + R_{1} R_{2} s^{2} Q (0) + R_{1} R_{3} s^{2} Q (0) + R_{2} R_{3} s^{2} Q (0) + L R_{2} s Q (0) + L R_{3} s Q (0))}{(s^{2} + 1) (R_{1} + R_{2} + L s + C L R_{2} s^{2} + C L R_{3} s^{2} + C R_{1} R_{2} s + C R_{1} R_{3} s + C R_{2} R_{3} s)}$

I1_LT와 Q_LT의 라플라스 역변환을 계산하여 $I_{1}$ 과 $Q$ 를 계산합니다. 결과를 단순화합니다. 출력값이 너무 길기 때문에 이를 표시하지 않습니다.

I1sol = ilaplace(I1_LT,s,t);
Qsol = ilaplace(Q_LT,s,t);
I1sol = simplify(I1sol);
Qsol = simplify(Qsol);

값 대입하기

결과를 플로팅하기 전에 먼저 기호 변수에 회로 소자의 수치적 값을 대입합니다. $R_{1} = 4 Ω$ , $R_{2} = 2 Ω$ , $R_{3} = 3 Ω$ , $C = 1 / 4 F$ , $L = 1.6 H$ 라고 하겠습니다. 초기 전류는 $I_{1} (0) = 2 A$ 이고 초기 전하는 $Q (0) = 2 C$ 라고 가정합니다.

vars = [R L C I1(0) Q(0)];
values = [4 2 3 1.6 1/4 2 2];
I1sol = subs(I1sol,vars,values)

I1sol = 
 $\frac{200 \cos (t)}{8161} + \frac{405 \sin (t)}{8161} + \frac{16122 e^{- \frac{81 t}{40}} (\cosh (\frac{\sqrt{1761} t}{40}) - \frac{742529 \sqrt{1761} \sinh (\frac{\sqrt{1761} t}{40})}{14195421})}{8161}$

Qsol = subs(Qsol,vars,values)

Qsol = 
 $\frac{924 \sin (t)}{8161} - \frac{1055 \cos (t)}{8161} + \frac{17377 e^{- \frac{81 t}{40}} (\cosh (\frac{\sqrt{1761} t}{40}) + \frac{1109425 \sqrt{1761} \sinh (\frac{\sqrt{1761} t}{40})}{30600897})}{8161}$

결과 플로팅하기

전류 I1sol과 전하 Qsol을 플로팅합니다. 2개의 서로 다른 시간 구간 $0 \leq t \leq 15$ 와 $2 \leq t \leq 25$ 를 사용하여 과도 상태 동작과 정상 상태 동작을 모두 표시합니다.

subplot(2,2,1)
fplot(I1sol,[0 15])                      
title('Current')
ylabel('I1(t)')
xlabel('t')

subplot(2,2,2)
fplot(Qsol,[0 15])                        
title('Charge')
ylabel('Q(t)')
xlabel('t')

subplot(2,2,3)
fplot(I1sol,[2 25])                  
title('Current')
ylabel('I1(t)')
xlabel('t')
text(3,-0.1,'Transient')
text(15,-0.07,'Steady State')

subplot(2,2,4)
fplot(Qsol,[2 25])                        
title('Charge')
ylabel('Q(t)')
xlabel('t')
text(3,0.35,'Transient')
text(15,0.22,'Steady State')

결과 분석하기

초기에는 전류와 전하가 지수적으로 감소합니다. 그러나 장기적으로는 진동합니다. 이러한 동작은 각각 "과도 상태"와 "정상 상태"라고 합니다. 기호 결과를 사용하면 결과의 속성을 분석할 수 있습니다. 이는 수치 결과를 사용할 때는 가능하지 않습니다.

I1sol과 Qsol을 시각적으로 검토합니다. 항들의 합을 볼 수 있습니다. children을 사용하여 항을 구합니다. 그런 다음 [0 15]에 대해 항을 플로팅하여 항의 비중을 구합니다. 플롯에서 과도 상태 항과 정상 상태 항을 볼 수 있습니다.

I1terms = children(I1sol);
I1terms = [I1terms{:}];
Qterms = children(Qsol);
Qterms = [Qterms{:}];

figure;
subplot(1,2,1)
fplot(I1terms,[0 15])
ylim([-0.5 2.5])
title('Current terms')

subplot(1,2,2)
fplot(Qterms,[0 15])
ylim([-0.5 2.5])
title('Charge terms')

Figure contains 2 axes objects. Axes object 1 with title Current terms contains 3 objects of type functionline. Axes object 2 with title Charge terms contains 3 objects of type functionline.

플롯을 통해 I1sol에는 하나의 과도 상태 항과 하나의 정상 상태 항이 있고 Qsol에는 하나의 과도 상태 항과 두 개의 정상 상태 항이 있음을 볼 수 있습니다. 시각적 검토를 통해 I1sol과 Qsol에 exp 함수를 포함하는 하나의 항이 있는 것을 알 수 있습니다. 이 항이 과도 지수 감쇠의 원인이 된다고 가정합니다. has를 사용하여 exp가 있는 항이 있는지 확인함으로써 I1sol과 Qsol의 과도 상태 항과 정상 상태 항을 분리합니다.

I1transient = I1terms(has(I1terms,'exp'))

I1transient = 
 $\frac{16122 e^{- \frac{81 t}{40}} (\cosh (\frac{\sqrt{1761} t}{40}) - \frac{742529 \sqrt{1761} \sinh (\frac{\sqrt{1761} t}{40})}{14195421})}{8161}$

I1steadystate = I1terms(~has(I1terms,'exp'))

I1steadystate = 
 $(\begin{array}{c} \frac{200 \cos (t)}{8161} & \frac{405 \sin (t)}{8161} \end{array})$

마찬가지로, Qsol을 과도 상태 항과 정상 상태 항으로 분리합니다. 이 결과는 기호 계산이 문제를 분석하는 데 어떻게 도움이 되는지 보여줍니다.

Qtransient = Qterms(has(Qterms,'exp'))

Qtransient = 
 $\frac{17377 e^{- \frac{81 t}{40}} (\cosh (\frac{\sqrt{1761} t}{40}) + \frac{1109425 \sqrt{1761} \sinh (\frac{\sqrt{1761} t}{40})}{30600897})}{8161}$

Qsteadystate = Qterms(~has(Qterms,'exp'))

Qsteadystate = 
 $(\begin{array}{c} - \frac{1055 \cos (t)}{8161} & \frac{924 \sin (t)}{8161} \end{array})$