RC circuit DE solution

[Federico Manfrin]

Hi, I'm writing a notebook to reproduce what I can see in this wikipedia page:

https://it.wikipedia.org/wiki/Circuito_RC

The problem is the CAP 4, I'm looking forward to get the result that follow that text:
dove K è una costante. Dunque: (latex folleow)

{\displaystyle v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}

This is my code, running in a colab notebook:

# define the independent variable ‘t’

from sympy.abc import t

import sympy as sp

import math

w = 2 * math.pi * 50

# define dependent variable in symbol form

vc = sp.symbols('vc', cls=sp.Function)

C, R, vo = sp.symbols('C R vo')

cos = sp.cos

diffeq = vc(t).diff(t) + (1/(R*C))*vc(t) - (vo*cos(w*t)/(R*C))

res = sp.dsolve(diffeq, ics={vc(0): 0})

res

The result is really confused, compared to the nice wikipedia result.

Is there a way to get the same result? (follow the latex)

{\displaystyle v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}

Thanks so much

[Aaron Meurer]

What's the original ODE you are trying to solve? The solution you give
has symbols in it that aren't in your diffeq.

One change I would recommend is to use sympy.pi instead of math.pi.
math.pi is a float approximation to pi, whereas sympy.pi is pi
exactly. Alternatively, you can use w = symbols('omega') if you want
the solution to match the general text version.

Aaron Meurer

> --
> You received this message because you are subscribed to the Google Groups "sympy" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to sympy+un...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/sympy/8f5df611-4aa2-4f09-9f5b-a9d528c24333n%40googlegroups.com.

[Federico Manfrin]

Hi Aaron, thank you for the reply.

The code I have in my notebook (I'm using google colab) is this now (at the end of the post), but using w is not working anymore.
Anyway, if I move back to use the line w = 2 * sp.pi * 50 I still get a unreadble result.

The result should be something like  the nice wikipedia result I found there:
https://it.wikipedia.org/wiki/Circuito_RC#Risposta_in_frequenza_del_circuito_RC

Is there a way to get the a readble result? (something like the following latex)

{\displaystyle v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}

The wikipedia result has tau=R*C, theta is an angle displacement,  and K don't know what is.
In the italian wikipedia there's this piece (I try to translate) that I would like to reproduce in a notebook:

Let's see how does the RC circuit works with a sine wave. We can use voltage Kirchhoff law:
{\displaystyle V_{0}\cos(\omega t)=R\cdot i(t)+v_{C}(t)}

we can rewrite the equation like:
{\displaystyle V_{0}\cos(\omega t)=RC{\frac {dv_{C}(t)}{dt}}+v_{C}(t)}

and then solve the differential equation with constant coefficients with a known therm:
{\displaystyle {\frac {dv_{C}(t)}{dt}}+{\frac {1}{\tau }}v_{C}(t)={\frac {V_{0}\cos(\omega t)}{\tau }}}

where \tau =RC is still the time constant of the circuit.
The general solution come from the sum of the associated homogeneous solution:

v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}

and a particolar solution:
{\displaystyle K\cos(\omega t+\theta )\ }

where K is a constant. So:

{\displaystyle v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}

---------------- MY CODE --------------------

# define the independent variable ‘t’

from sympy.abc import t

import sympy as sp

#w = 2 * sp.pi * 50

w = sp.symbols('omega')

# define dependent variable in symbol form

vc = sp.symbols('vc', cls=sp.Function)

C, R, vo = sp.symbols('C R vo')

cos = sp.cos

diffeq = vc(t).diff(t) + (1/(R*C))*vc(t) - (vo*cos(w*t)/(R*C))

res = sp.dsolve(diffeq, ics={vc(0): 0})

res

Il giorno lunedì 27 giugno 2022 alle 22:30:32 UTC+2 ... ha scritto:

What's the original ODE you are trying to solve? The solution you give
has symbols in it that aren't in your diffeq.

One change I would recommend is to use sympy.pi instead of math.pi.
math.pi is a float approximation to pi, whereas sympy.pi is pi
exactly. Alternatively, you can use w = symbols('omega') if you want
the solution to match the general text version.

Aaron Meurer

On Mon, Jun 27, 2022 at 2:14 PM FM

[Oscar]

On Tuesday, 28 June 2022 at 08:03:05 UTC+1 federic...@gmail.com wrote:

Hi Aaron, thank you for the reply.

The code I have in my notebook (I'm using google colab) is this now (at the end of the post), but using w is not working anymore.
Anyway, if I move back to use the line w = 2 * sp.pi * 50 I still get a unreadble result.

The result should be something like  the nice wikipedia result I found there:
https://it.wikipedia.org/wiki/Circuito_RC#Risposta_in_frequenza_del_circuito_RC

Is there a way to get the a readble result? (something like the following latex)

{\displaystyle v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}

The wikipedia result has tau=R*C, theta is an angle displacement,  and K don't know what is.

The result that I see is similar:

In [13]: diffeq

Out[13]:
d           vo⋅cos(100⋅π⋅t)   vc(t)
──(vc(t)) - ─────────────── + ─────

dt                C⋅R          C⋅R

In [12]: res
Out[12]:
                                                                  -t      
                                                                  ───      
                                                                  C⋅R      
        100⋅π⋅C⋅R⋅vo⋅sin(100⋅π⋅t)    vo⋅cos(100⋅π⋅t)          vo⋅ℯ        
vc(t) = ───────────────────────── + ────────────────── - ──────────────────
                   2  2  2                 2  2  2              2  2  2    
            10000⋅π ⋅C ⋅R  + 1      10000⋅π ⋅C ⋅R  + 1   10000⋅π ⋅C ⋅R  + 1

The main difference is that the constant K is given explicitly here in terms of C, R and v0. If you want to substitute for tau and K you can do something like:

In [19]: tau, K = symbols('tau, K')

In [20]: res.subs(C, tau/R).subs(vo, K*(1 + 10000*pi**2*tau**2))
Out[20]:
                                                     -t
                                                     ───
                                                      τ
vc(t) = 100⋅π⋅K⋅τ⋅sin(100⋅π⋅t) + K⋅cos(100⋅π⋅t) - K⋅ℯ  

You can see now something more similar. The difference is that you have sin and cos rather than cos with a phase shift theta but those are equivalent under trigonometric identities. Another difference is the constant K on the exponential term but the wikipedia solution is incorrect there:

 

In the italian wikipedia there's this piece (I try to translate) that I would like to reproduce in a notebook:

Let's see how does the RC circuit works with a sine wave. We can use voltage Kirchhoff law:
{\displaystyle V_{0}\cos(\omega t)=R\cdot i(t)+v_{C}(t)}

we can rewrite the equation like:
{\displaystyle V_{0}\cos(\omega t)=RC{\frac {dv_{C}(t)}{dt}}+v_{C}(t)}

and then solve the differential equation with constant coefficients with a known therm:
{\displaystyle {\frac {dv_{C}(t)}{dt}}+{\frac {1}{\tau }}v_{C}(t)={\frac {V_{0}\cos(\omega t)}{\tau }}}

where \tau =RC is still the time constant of the circuit.
The general solution come from the sum of the associated homogeneous solution:

v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}

It is incorrect to use v_C(0) here. There should be a constant there but it is not generally equal to the initial condition.

and a particolar solution:
{\displaystyle K\cos(\omega t+\theta )\ }

where K is a constant. So:

{\displaystyle v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}

You can see that this is incorrect if you just substitute t=0.

 

--

Oscar

[Federico Manfrin]

Thank you Oscar,

onestly I don't get the same resoult of yours, because I get 3 result  , not just one like you:

one solution for for τ/R ​=−i/(100πRi) , one solution for τ/R ​=i/(100πRi) and another is otherwise

Anyway, the otherwise seems to be possible to be simplified and become like your solution, so thank you very much for the help.

It should be nice to avoid the 3 solution and keep just the otherwise one and get it simplified, but I just want to thank you for your help.

[Michael Gustafson]

I know this is almost a year old, but is this notebook helpful to you:

https://colab.research.google.com/drive/1p5WP-twoxD5HKq9gziNWKqSi8Xo_cotg?usp=sharing

-Michael

[Michael Gustafson]

Also, here is the code - apologies for not including it before:

This notebook was prompted by the thread at:<br>
https://groups.google.com/g/sympy/c/VW2eFNBw1Cg/m/AxTlTwGUHQAJ


import sympy as sp

sp.init_printing(use_latex=True)

'''
Without the real=True part, there are some solutions with imaginary components
'''

w, t, tau, vo, vc0, K, theta = sp.var(r'\omega, t, \tau, v_o, v_{C}(0), K, \theta', real=True)

vc = sp.Function('v_c')

diffeq = sp.Derivative(vc(t), t) + vc(t) / tau  - vo * sp.cos(w * t) / tau
diffeq

res = sp.dsolve([diffeq], [vc(t)], ics={vc(0):vc0})
res

'''

Quick note - the solution at<br>
https://it.wikipedia.org/wiki/Circuito_RC#Risposta_in_frequenza_del_circuito_RC
<br>is not quite right in that<br>
$$ {\displaystyle v_{C}(t)=v_{C}(0)e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}$$
gets the initial condition wrong; I think it should be something like:
$${\displaystyle v_{C}(t)=(v_{C}(0)-K\cos(\theta))e^{-{\frac {t}{\tau }}}+K\cos(\omega t+\theta )}$$
'''

# Collecting vo terms
res[0].rhs.expand().collect(vo).trigsimp()

# Collecting exponential terms
res[0].rhs.expand().collect(sp.exp(-t/tau)).trigsimp()

'''
The $K$ and $\theta$ here match https://it.wikipedia.org/wiki/Circuito_RC#Metodo_simbolico_per_la_risposta_in_frequenza (using trig identities to get single magnitude and phase)
'''

K = vo * sp.sqrt((w*tau)**2 + 1**2) / ((w*tau)**2 + 1)
K

theta = sp.atan2(-w*tau, 1)
theta

vcneat = (vc0-K*sp.cos(theta))*sp.exp(-t/tau)+K*sp.cos(w*t+theta)
vcneat

# Are they the same?
(vcneat - res[0].rhs).trigsimp()

[Federico Manfrin]

Hi Michael,
I'm back in this notebook for teaching purpose, and I really want to thank you for the notebook you shared.

Thank you
Federico