This notebook simulates the pulse response of a basic RLC circuit. You can use this notebook to help you choose a resistor R that you solder in series at your pulse-generating source, to minimize ringing overshoot (due to inductance) while not rounding the pulse edges too much (due to capacitance). In the last plot below, we find a good resistance value for R.
Define the function that performs an inverse Laplace transform¶
# PULSE_RESPONSE_EXAMINER PROGRAM # Version 1.0.1, 02-Jul-2021 # Copyright © 2021 Joe Czapski # Contact: xejgyczlionoldmailservertigernet replace lion with # an at sign and tiger with a dot. # # Licensed under the Apache License, Version 2.0 (the "License"); # you may not use this file except in compliance with the License. # You may obtain a copy of the License at # # http://www.apache.org/licenses/LICENSE-2.0 # # Unless required by applicable law or agreed to in writing, software # distributed under the License is distributed on an "AS IS" BASIS, # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. # See the License for the specific language governing permissions and # limitations under the License. import math, cmath import numpy as np from numpy.fft import fft import matplotlib.pyplot as plt # FUNCTION inv_laplace(): employs the 'fourier method' for inversion of the Laplace # transform. The algorithm here is an adaptation by Joe Czapski of a Mathcad # worksheet by Sven Lindhardt. It's a practical approximation of inverse Laplace for # systems close to stability, whether a little overstable or understable, when the # rightmost pole is not greater than zero. NumPy's complex FFT is used in this # application to generate a fast inverse Laplace transform. def inv_laplace(n, R, L, C, t1, t2, tspan): # Error correction factor set to 2.5 to reduce artifacts when both L and C are tiny. sigma = 2.5 / tspan # Calculate delta omega (angular frequency sample interval). domega = math.pi / tspan # Convert impedance and time values to complex numbers for transfer function. Rc = R+0.0j; Lc = L+0.0j; Cc = C+0.0j; t1c = t1+0.0j; t2c = t2+0.0j Vouts = np.zeros(n+1, dtype=np.cdouble) for i in range(0, n+1): omega = i * domega s = sigma - omega*1.0j sL = s * Lc sCinv = (1.0+0.0j) / (s * Cc) # Calculate H(s) for the desired circuit. # For an RLC circuit H(s) = (1/sC) / (R + sL + 1/sC) Hs = sCinv / (Rc + sL + sCinv) # Calculate Vin(s) for the input waveform. # For an ideal pulse Vin(s) = 100 * (e^-t1s - e^-t2s) / s Vins = (100.0+0.0j) * (cmath.exp(-t1c*s) - cmath.exp(-t2c*s)) / s # Calculate Vout(s) = H(s) * Vin(s) Vouts[i] = Hs * Vins # Perform the complex fourier transform. Before calling fft(), # construct a symmetric version of complex array Vouts. Vouts[n] = Vouts[n].real + 0.0j Vouts_symmetric = np.concatenate([Vouts, (Vouts[1:n])[::-1]]) Vout_unscaled = fft(Vouts_symmetric) # Take the real part of Vout_unscaled, first half of the symmetric array, # and scale for inverse Laplace using time span, sigma, and exp function. # The sample interval in seconds = tspan / n sigdx = sigma * tspan / n Vout = np.zeros(n) for i in range(0, n): Vout[i] = math.exp(i*sigdx) * Vout_unscaled[i].real / tspan return Vout
R = 10.0 # ohms L = 4.0e-9 # henries C = 5.0e-12 # farads t1 = 1.0e-8 # pulse rising edge seconds from t0 t2 = 2.0e-8 # pulse falling edge seconds from t0 # Set tspan which is the time span of the analysis in seconds from t0. # Artifacts can arise when Vout doesn't settle to zero before the end of the # analysis time. You may need to adjust tspan to increase the analysis time # following the end of the pulse to allow settling to zero. tspan = 4.0e-8 # Set n which is the number of points to use for the analysis. n must be a power of 2. # Set n=2048 for most cases. Set n=8192 to reduce artifacts when both L and C are tiny. n = 2048 ns_per_s = 1.0e9 # time display unit correction xdata = np.linspace(0.0, tspan * ns_per_s * (n-1)/n, n) ydata = inv_laplace(n, R, L, C, t1, t2, tspan) plt.rcParams['figure.figsize'] = [12, 6] plt.title('Output Voltage Waveform, R = %d Ω' % (R)) plt.xlabel('Time ns') plt.ylabel('Magnitude %') plt.grid(True) plt.plot(xdata, ydata, color='firebrick') plt.show()
R = 200.0 # ohms ydata = inv_laplace(n, R, L, C, t1, t2, tspan) plt.title('Output Voltage Waveform, R = %d Ω' % (R)) plt.xlabel('Time ns') plt.ylabel('Magnitude %') plt.grid(True) plt.plot(xdata, ydata, color='firebrick') plt.show()
R = 50.0 # ohms ydata = inv_laplace(n, R, L, C, t1, t2, tspan) plt.title('Output Voltage Waveform, R = %d Ω' % (R)) plt.xlabel('Time ns') plt.ylabel('Magnitude %') plt.grid(True) plt.plot(xdata, ydata, color='firebrick') plt.show()