diff sim.py @ 10:2832aefd442c

Add 'no driver' version which can simulate in normal mode. Add Python code to optimise design.
author Daniel O'Connor <darius@dons.net.au>
date Fri, 17 Nov 2023 23:47:05 +1030
parents
children 8d34a9eec184
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/sim.py	Fri Nov 17 23:47:05 2023 +1030
@@ -0,0 +1,148 @@
+import datetime
+import pathlib
+from scipy.optimize import differential_evolution
+from spicelib.simulators.ltspice_simulator import LTspice
+from spicelib.log.ltsteps import LTSpiceLogReader
+from spicelib.editor.spice_editor import SpiceEditor
+#from spicelib.sim.sim_runner import SimRunner
+#import shlex
+import sqlite3
+import subprocess
+
+def simulate(src, simexe, simflags, rundir, runname, params):
+    src = pathlib.Path(src)
+    runname = pathlib.Path(runname)
+    assert(runname.suffix == '.net')
+    rundir = pathlib.Path(rundir)
+
+    e = SpiceEditor(src)
+    for k in params:
+        e.set_parameter(k, params[k])
+
+    #s = SimRunner(simulator = LTspice, output_folder = 'tmp', verbose = True)
+    #s.simulator.spice_exe = [simexe]
+    #raw_file, log_file = s.run_now(e, run_filename = runname)
+
+    e.write_netlist(rundir / runname)
+
+    cmd = [simexe, '-b']
+    cmd.extend(simflags)
+    cmd.append(runname.name)
+    #print(' '.join(map(shlex.quote, cmd)))
+    then = datetime.datetime.now()
+    p = subprocess.Popen(cmd, cwd = rundir)
+    p.communicate()
+    now = datetime.datetime.now()
+    taken = now - then
+
+    if p.returncode != 0:
+        raise Exception('failed')
+
+    #rawpath = rundir / runname.with_suffix('.raw')
+    #raw = spicelib.RawRead(rawpath)
+    logpath = rundir / runname.with_suffix('.log')
+    log = LTSpiceLogReader(logpath)
+
+    power = log.get_measure_value('pout')
+    eff = log.get_measure_value('efficiency')
+    thd = log.fourier['V(rfout)'][0].thd
+    ipeak_u2 = log.get_measure_value('ipeak_u2')
+    ipeak_u5 = log.get_measure_value('ipeak_u5')
+    ipeak = max(ipeak_u2, ipeak_u5)
+
+    return then, taken, power, eff, thd, ipeak
+
+
+def fn(v, cur, simexe, simflags, rundir, circ, ctr):
+    '''Called by differential_evolution, v contains the evolved parameters, the rest are passed in from the args parameter'.
+    Returns the cost value which is being minimised'''
+
+    # Get parameters for run
+    duty, c1, c2, l1, l2 = v
+
+    cur.execute('SELECT cost, run, power, efficiency, thd, ipeak FROM GAN190 WHERE duty = ? AND c1 = ? AND c2 = ? AND l1 = ? AND l2 = ?', (duty, c1, c2, l1, l2))
+    tmp = cur.fetchone()
+    if tmp is not None:
+        cost, run, power, eff, thd, ipeak = tmp
+        print(f'Found run {run:3d}: Duty {duty:3.0f}%, C1 {c1:3.0f}pF, C2 {c2:3.0f}pF, L1 {l1:3.0f}uH, L2 {l2:4.0f}nH -> Power: {power:6.1f}W Efficiency: {eff:5.1f}% THD: {thd:5.1f}% IPeak: {ipeak:4.1f}A Cost: {cost:6.2f}')
+        return cost
+
+    # Get next run number
+    run = next(ctr)
+
+    # Run the simulation
+    # Need to convert units to suit
+    runname = pathlib.Path(circ).stem + '-' + str(run) + '.net'
+    then, taken, power, eff, thd, ipeak = simulate(circ, simexe, simflags, rundir, runname,
+                                                   {'dutypct' : duty, 'C1' : c1 * 1e-12, 'C2' : c2 * 1e-12, 'L1' : l1 * 1e-6, 'L2' : l2 * 1e-9})
+
+    # Calculate the cost
+    tpwr = 2000
+    tthd = 2
+    teff = 90
+    imax = 11
+
+    cost = 0
+    if power < tpwr / 2:
+        cost += 100
+    elif power < tpwr:
+        cost += (tpwr - power) / tpwr / 100
+
+    if thd > 5 * tthd:
+        cost += 100
+    else:
+        thdinv = 100 - thd
+        tthdinv = 100 - tthd
+        if thdinv < tthdinv:
+            cost += (tthdinv - thdinv) / tthdinv
+    if eff < teff:
+        cost += (teff - eff) / teff
+
+    if ipeak > imax:
+        cost += (ipeak - imax) * 2
+
+    # Log & save the results
+    print(f'Run {run:3d}: Duty {duty:3.0f}%, C1 {c1:3.0f}pF, C2 {c2:3.0f}pF, L1 {l1:3.0f}uH, L2 {l2:4.0f}nH -> Power: {power:6.1f}W Efficiency: {eff:5.1f}% THD: {thd:5.1f}% IPeak: {ipeak:4.1f}A Cost: {cost:6.2f}')
+    taken = taken.seconds + taken.microseconds / 1e6
+    cur.execute('INSERT INTO GAN190 (name, run, start, duration, duty, c1, c2, l1, l2, power, efficiency, thd, ipeak, cost) VALUES (?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?, ?)',
+                 (runname, run, then, taken, duty, c1, c2, l1, l2, power, eff, thd, ipeak, cost))
+    cur.execute('COMMIT')
+
+    return cost
+
+def ev():
+    # Bounds for parameters
+    bounds = [(10, 80),
+              (1, 20),
+              (10, 300),
+              (1, 20),
+              (10, 500)]
+    # Initial solution
+    x0 = [36, 10, 155, 5, 140]
+    dbh = sqlite3.connect('results.db')
+    cur = dbh.cursor()
+    cur.execute('''
+CREATE TABLE IF NOT EXISTS GAN190 (
+    name	TEXT,		-- Circuit name
+    run		INTEGER,	-- Run number
+    start	DATETIME,	-- Datetime run started
+    duration	REAL,		-- Length of run (seconds)
+    duty	REAL,		-- Duty cyle (%)
+    c1		REAL,		-- Value of C1 (pF)
+    c2		REAL,		-- Value of C2 (pF)
+    l1		REAL,		-- Value of L1 (uH)
+    l2		REAL,		-- Value of L2 (nH)
+    power	REAL,		-- Measured power (Watts)
+    efficiency	REAL,		-- Measured efficiency (%)
+    thd		REAL,		-- Total harmonic distortion (%)
+    ipeak	REAL,		-- Peak drain current (A)
+    cost	REAL		-- Calculated cost metric
+);''')
+    def ctr():
+        i = 0
+        while True:
+            yield i
+            i += 1
+    return differential_evolution(fn, bounds, x0 = x0,
+                                  args = (cur, '/Users/oconnd1/bin/runltspice', ['-alt'], 'tmp', 'pa-GAN190-PP-nodriver.net', ctr()),
+                                  integrality = True, disp = True, seed = 12345)