Modelling a diode for use in SPICE simulations

I have recently been playing around with LTSPICE to gain a better understanding of some of the circuits I have been building. I got involved in looking a the I-V characteristics of some LEDs that I own and was trying to model them accurately. This led me into trying to produce my own SPICE parameters. I found a page that discussed the basics and had a downloadable spreadsheet that utilised an iterative approach to modelling diodes. This gave me an idea to utilise Python to perform a non-linear least squares regression of diode data to the model presented here. I have gone for the full explicit solution for the forward bias characteristics of a diode. The program below imports a set of I-V data (this can be generated from a datasheet using the WebPlotDigitizer chrome app) and fits it to a three parameter model to yield values for IS, N and VS (see http://www.allaboutcircuits.com/vol_3/chpt_3/14.html for an explanation of the parameters).

Here’s the program, it requires Python3 with the Numpy, Scipy, and Matplotlib modules. It can be run directly from the command line and just requires the name of the file containing the I-V data. There are a number of options that can be accessed via the –help flag. I hope someone finds it useful.

#!/usr/bin/env python3
import numpy as np
import matplotlib.pyplot as plt
import argparse
from scipy.optimize import curve_fit
from scipy.special import lambertw

def func(VS, IS, N, RS):
#Define model for diode (see wikipedia article)
    VT = 26e-3
    w = lambertw((IS * RS /(N * VT)) * np.exp((VS + IS * RS) / (N * VT)))
    Current = np.log10((IS * ((w * N * VT / (RS * IS)) - 1) * 1000))
    return Current.real

def plot(xdata, ydata, nPoints, IS, N, RS):
#generate a points to plot
    vMin = xdata[0]
    vMax = xdata[len(xdata) - 1],
    vRange = vMax - vMin
    vMin = vMin - 0.1 * vRange
    vMax = vMax + 0.1 * vRange
    vStep = (vMax - vMin) / (nPoints - 1)
    VS = []
    ID = []
    for i in range(0, nPoints):
        voltage = vMin + i * vStep
        current = func(voltage, IS, N, RS)
        VS.append(voltage)
        ID.append(10**current)

#Plot the data and model
    plt.figure().canvas.set_window_title('Plot Window') 
    plt.semilogy(xdata, ydata, 'r*', VS, ID, 'b-')
    plt.ylabel('Current / mA')
    plt.xlabel('Voltage / V')
    plt.title('Diode I-V Characteristic')
    plt.show()

def main():
#set up parser for command line args
    parser = argparse.ArgumentParser(prog='DiodeModel')
    parser.add_argument('filename', help='Name of file containing I-V data (I in mA, V in volts)')
    parser.add_argument('-c', '--convert', help='Convert read in current to mA', action="store_true")
    parser.add_argument('-p', '--plot', help='Just plot the data and initial guess, no fitting performed', action="store_true")
    parser.add_argument('-IS', '--IS', type=float, default=1e-14, help='Initial guess at saturation current (default = 1e-14 A')
    parser.add_argument('-N', '--N', type=float, default=1, help='Initial guess at Emission coefficient (default = 1)')
    parser.add_argument('-RS', '--RS', type=float, default=10, help='Initial guess at ohmic resistance (default = 10 ohm)')
    parser.add_argument('-m', '--maxit', type=int, default=1000, help='Maximum number of iterations (default = 1000)')
    parser.add_argument('-n', '--npoints', type=int, default=250, help='Number of points in plot (default = 250)')
    parser.add_argument('-s', '--save', help='Save fit parameters to file', action="store_true")
    parser.add_argument('-f', '--fitfile', type=str, default="", nargs='?', help='Load fit parameters from specified file')
    args = parser.parse_args()
    
#Read in data from file (V in volt, I in milliamp)
    xdata, ydata = np.loadtxt(args.filename, unpack=True)
    if args.convert:           #Data is in A convert to mA
        ydata = ydata * 1000
    logydata = np.log10(ydata) #log of current to produce stable fit

#Set up initial guess
    params = dict(IS = args.IS, N = args.N, RS = args.RS)

#Check if input file is given
    if args.fitfile is None:                              #Use default filename
        fitFile = args.filename.split('.')[0] + ".fit"
        try:
            fh = open(fitFile, 'r')
            for line in fh:
                name = line.split('=')[0].strip()
                value = float(line.split('=')[1].strip())
                params[name] = value
            fh.close()
        except IOError as e:
            print("I/O error({0}): {1}".format(e.errno, e.strerror))
    elif len(args.fitfile) > 0:                           #Use input filename
        fitFile = args.fitfile
        try:
            fh = open(fitFile, 'r')
            for line in fh:
                name = line.split('=')[0].strip()
                value = float(line.split('=')[1].strip())
                params[name] = value
            fh.close()
        except IOError as e:
            print("I/O error({0}): {1}".format(e.errno, e.strerror))
    
    if args.plot:
        try:
#Plot data and initial guess
            print("Plotting characteristic with following parameters" + 
                "\nIS = " + str(params['IS']) + "\nN = " + str(params['N']) + "\nRS = " + str(params['RS']))
            plot(xdata, ydata, args.npoints, params['IS'], params['N'], params['RS'])
        except:
            print("Plotting error")
            exit()
    else:
#Perform non-linear least squares fit
        try:
            popt, pcov, infodict, errmsg, ier = curve_fit(
                func, xdata, logydata, p0=(params['IS'], params['N'], params['RS']), maxfev = args.maxit, full_output = True)
            
#Print converged fit parameters
            print("Fit converged in " + str(infodict['nfev']) + " iterations with the following parameters")
            printString = "IS = " + str(popt[0]) + "\nN = " + str(popt[1]) + "\nRS = " + str(popt[2])
            print(printString)

#Output parameters to file
            if args.save:
                print("\nWriting fit parameters to file")
                outFile = args.filename.split('.')[0] + ".fit"
                try:
                    fh = open(outFile,"w")
                    fh.write(printString)
                    fh.close()
                except IOError as e:
                    print("I/O error({0}): {1}".format(e.errno, e.strerror))

#generate a plot of the fit
            plot(xdata, ydata, args.npoints, popt[0], popt[1], popt[2])

		except RuntimeError as e:
			print("Error ā€“ Fit did not converge, try adjusting the starting guess")
			print(e)

if __name__ == "__main__":
    main()
Advertisements

13 thoughts on “Modelling a diode for use in SPICE simulations

  1. I’d love to see an example input file and program invocation. I am new to python.

    Otherwise, this looks like a cool piece of software.

    • Hi, just open Notepad or similar and type your values like this

      1.0 100.01
      1.5 200.46
      2.0 300.78

      then save this as txt or whatever you want and provide this values to the script.

      Example
      python scriptname.py filename.txt

  2. Excellent script and clever use of WebPlotDigitizer!

    I followed your instructions with the IV graph given on the stack exchange page you linked: http://electronics.stackexchange.com/questions/9510/how-do-i-model-an-led-with-spice

    … and WebPlotDigitizer gave me the following values (had to remove the commas from the generated data.csv file):

    1.47159725771882 0.1675368758897402
    1.5417367942905789 0.44016347920521
    1.5583259886891418 0.7562060928264692
    1.599611421691553 1.4727542814986494
    1.624408040771061 2.4106130614445678
    1.6986790263936722 6.8422380316915845
    1.847083634882728 33.34501053719433
    1.9459055715786266 63.04666099093824
    2.048779709581435 96.36361760157867
    2.260637986644994 183.74845113160293

    However, and with the default values for IS,N,RS, your script came up with a negative IS value.

    Adding this test as the first line of your func(VS, IS, N, RS) solved the problem:
    if IS < 0 or N < 0 or RS < 0: return 1e10

    With my bodge, these are the values I get:
    Fit converged in 249 iterations with the following parameters
    IS = 7.6000773057e-16
    N = 2.17306361372
    RS = 2.14797836067

    The fit is smooth with no visible cusp. Also tested with Python 2.7.9, your script works fine. Hope this helps šŸ™‚

    Cheers,
    John

  3. Great Job ! In combination with WebPlotDigitizer, Chrome and Spice the usefuliest App I saw for years. You should correct your code above, because Python won’t accept HTML stuff like &qout;

  4. Could you put the code on github ?
    As is, your code can’t be used directly because the HTML formatting ruins
    a whole bunch of special chars

  5. Nice! I have found one little problem though, if maxfev is reached, curve_fit raises a RuntimeError therefore the variable ‘errmsg’ doesn’t have any value, causing an UnboundLocalError at line 123. I solved in a simple way, on lines 121-123

    except RuntimeError as e:
    print(“Error – Fit did not converge, try adjusting the starting guess”)
    print(e)

    In addition WordPress put a > on line 75.
    Cheers,

    Davide

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s