Skip to content

Instantly share code, notes, and snippets.

@rezamarzban

rezamarzban/RLC.py

Last active October 26, 2025 05:38
Show Gist options

  • Select an option

    Learn more about clone URLs
  • Save rezamarzban/c8a392427357116bd68034d6763c8f5d to your computer and use it in GitHub Desktop.

Select an option

Learn more about clone URLs
Save rezamarzban/c8a392427357116bd68034d6763c8f5d to your computer and use it in GitHub Desktop.
Download ZIP
Raw
RLC.py
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
# Parameters
R = 0.01
L = 1e-6
C = 1e-7
f = 50e3 # 50 kHz
omega = 2 * np.pi * f
T_period = 1 / f # 20 µs period
# State-space function for RLC
def rlc(y, t):
i, vc = y
vin = 2.5 + 2.5 * np.sign(np.sin(omega * t)) # Square wave 0-5 V, starts at 2.5 + 2.5 = 5 V
didt = (vin - R * i - vc) / L
dvcdt = i / C
return [didt, dvcdt]
# Simulate full transient from t=0 to 10 ms (sufficient for settling)
t_full_start = 0
t_full_end = 10e-3
dt_full = 10e-9 # Fine steps
t_full = np.arange(t_full_start, t_full_end + dt_full, dt_full)
# Initial conditions: zero (uic equivalent)
y0 = [0, 0]
# Solve ODE for full time
sol = odeint(rlc, y0, t_full)
vc_full = sol[:, 1]
# Define windows: 200 µs segments at different starting times
window_duration = 200e-6
start_times = [0, 1e-3, 5e-3, 9e-3] # 0 µs, 1 ms, 5 ms, 9 ms
# Create subplots
fig, axes = plt.subplots(len(start_times), 1, figsize=(12, 4 * len(start_times)), sharex=True)
if len(start_times) == 1:
axes = [axes]
for idx, t_window_start in enumerate(start_times):
ax = axes[idx]
# Find indices for this window
idx_start = np.argmin(np.abs(t_full - t_window_start))
idx_end = np.argmin(np.abs(t_full - (t_window_start + window_duration)))
t_window = t_full[idx_start:idx_end + 1] - t_window_start # Relative time
vc_window = vc_full[idx_start:idx_end + 1]
# Compute v_in for window (for reference)
vin_window = 2.5 + 2.5 * np.sign(np.sin(omega * (t_window + t_window_start)))
# Plot
ax.plot(t_window * 1e6, vin_window, 'k--', linewidth=1, alpha=0.7, label='v_in(t)')
ax.plot(t_window * 1e6, vc_window, 'b-', linewidth=1, label='v_c(t)')
# Compute ripple in this window
ripple = np.ptp(vc_window)
ax.set_ylabel('Voltage (V)')
ax.set_title(f'Start Time: {t_window_start * 1e3:.0f} ms\nRipple: {ripple:.2f} V')
ax.grid(True, alpha=0.3)
ax.legend()
ax.set_ylim(-1, 8)
if idx == len(start_times) - 1:
ax.set_xlabel('Time in Window (µs)')
plt.tight_layout()
plt.show()
# Print summary
print("Summary of Ripple in Each 200 µs Window:")
print("-" * 50)
for t_window_start in start_times:
idx_start = np.argmin(np.abs(t_full - t_window_start))
idx_end = np.argmin(np.abs(t_full - (t_window_start + window_duration)))
ripple = np.ptp(vc_full[idx_start:idx_end + 1])
print(f"Start: {t_window_start * 1e3:.0f} ms → Ripple: {ripple:.2f} V")
Raw
RLC_tss.py
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
def calculate_steady_state_timing(R=0.01, L=1e-6, C=1e-7, f=50e3, epsilon=0.01):
"""
Calculate the steady-state timing for an RLC circuit.
Parameters:
- R: Resistance in ohms (default: 0.01)
- L: Inductance in henries (default: 1e-6)
- C: Capacitance in farads (default: 1e-7)
- f: Clock frequency in Hz (default: 50e3) - NOTE: Does not affect settling time
- epsilon: Tolerance level for transient decay (default: 0.01 for 1%)
Returns:
- Dictionary with alpha, tau, t_ss, zeta, is_underdamped
"""
alpha = R / (2 * L)
tau = 1 / alpha
t_ss = -np.log(epsilon) / alpha
zeta = R / (2 * np.sqrt(L / C))
is_underdamped = zeta < 1
return {
'alpha': alpha,
'tau': tau,
't_ss': t_ss,
'zeta': zeta,
'is_underdamped': is_underdamped
}
# State-space function for RLC
def rlc(y, t, R, L, C, omega):
i, vc = y
vin = 2.5 + 2.5 * np.sign(np.sin(omega * t)) # Square wave 0-5 V
didt = (vin - R * i - vc) / L
dvcdt = i / C
return [didt, dvcdt]
# Main function to plot response up to steady-state time
def plot_response_to_steady_state(R=0.01, L=1e-6, C=1e-7, f=50e3, epsilon=0.01):
# Calculate steady-state parameters
result = calculate_steady_state_timing(R, L, C, f, epsilon)
t_ss = result['t_ss']
alpha = result['alpha']
zeta = result['zeta']
print(f"Damping factor (alpha): {alpha:.2f} s^-1")
print(f"Damping ratio (zeta): {zeta:.5f}")
print(f"Steady-state time (t_ss): {t_ss * 1e3:.2f} ms")
# Time array: from 0 to 1.5 * t_ss for visualization
t_start = 0
t_end = 1.5 * t_ss
dt = 10e-9 # Fine resolution
t = np.arange(t_start, t_end + dt, dt)
omega = 2 * np.pi * f
# Initial conditions
y0 = [0, 0]
# Solve ODE
sol = odeint(rlc, y0, t, args=(R, L, C, omega))
vc = sol[:, 1]
# Compute v_in for plot
vin = 2.5 + 2.5 * np.sign(np.sin(omega * t))
# Plot
plt.figure(figsize=(12, 6))
plt.plot(t * 1e3, vin, 'k--', linewidth=1, alpha=0.7, label='v_in(t)')
plt.plot(t * 1e3, vc, 'b-', linewidth=1, label='v_c(t)')
plt.axvline(t_ss * 1e3, color='r', linestyle='--',
label=f'Steady-state at {t_ss * 1e3:.2f} ms (ε={epsilon})')
plt.xlabel('Time (ms)')
plt.ylabel('Voltage (V)')
plt.title('RLC Circuit Response Up to Steady-State Time')
plt.grid(True, alpha=0.3)
plt.legend()
plt.ylim(-1, 8)
plt.show()
# Example usage
plot_response_to_steady_state()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment