From 2f26f7709f846fcf08c503e200b922a5bc9a9053 Mon Sep 17 00:00:00 2001 From: Mattia Fiumara Date: Sun, 10 May 2026 14:48:34 +0200 Subject: [PATCH] feat(core): mutual inductance (K), per-element ic=, and .tran uic MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Adds the SPICE K-element for magnetically-coupled inductors, ic= initial conditions on capacitors and inductors, and the .tran ... uic flag that seeds the transient from those ICs instead of a DC operating point. Lets spice-ts run the Chua & Lin 8-7 (page 343) benchmark — a classic SPICE torture test where every simulator gives a slightly different answer because it stresses K-coupling, IC handling, and stiff multi-time-scale dynamics simultaneously. Wired into the showcase under "Benchmarks". UIC seed builds a substituted DC system (caps with ic become V sources, inductors with ic become I sources, K-elements drop) and maps the result back to the original circuit's branch layout. Inductor → I-source flips node order to match the inductor branch-current sign convention. The seeded solution does not satisfy the trapezoidal companion's history, so the driver bootstraps with Backward Euler on step 1 — same pattern used after waveform breakpoints. Tokenizer now ignores trailing unit letters after numeric values ("2H", "100Ohm", "1uF"); previously these threw ParseError. Issue #48. Co-Authored-By: Claude Opus 4.7 (1M context) --- examples/showcase/main.tsx | 32 +++++++ .../core/src/analysis/transient-driver.ts | 19 +++- packages/core/src/analysis/transient.ts | 2 + packages/core/src/analysis/uic.ts | 73 ++++++++++++++ packages/core/src/circuit.ts | 95 +++++++++++++++++-- packages/core/src/devices/capacitor.ts | 1 + packages/core/src/devices/index.ts | 1 + packages/core/src/devices/inductor.ts | 1 + .../core/src/devices/mutual-inductor.test.ts | 75 +++++++++++++++ packages/core/src/devices/mutual-inductor.ts | 43 +++++++++ packages/core/src/parser/index.ts | 42 ++++++-- packages/core/src/parser/tokenizer.ts | 14 ++- packages/core/src/simulate.ts | 7 +- packages/core/src/types.ts | 6 ++ 14 files changed, 391 insertions(+), 20 deletions(-) create mode 100644 packages/core/src/analysis/uic.ts create mode 100644 packages/core/src/devices/mutual-inductor.test.ts create mode 100644 packages/core/src/devices/mutual-inductor.ts diff --git a/examples/showcase/main.tsx b/examples/showcase/main.tsx index 50c0aef..ec1b4a1 100644 --- a/examples/showcase/main.tsx +++ b/examples/showcase/main.tsx @@ -205,6 +205,38 @@ D1 sw out DMOD C1 out 0 100u Rload out 0 10 .tran 50n 500u`, + }, + { + id: 'chua-lin', + name: 'Chua & Lin 8-7', + desc: 'Coupled inductors + UIC', + icon: '\u29c9', + group: 'Benchmarks', + tag: '.tran', + signals: ['x', 'y', 'out'], + tranNetlist: ` +* Chua & Lin "Computer-Aided Analysis" 8-7 (page 343). +* Three magnetically-coupled inductors (K1/K2/K3) with two capacitor and +* one inductor initial conditions. A classic SPICE benchmark \u2014 every +* simulator gives a slightly different answer because it forces the +* engine to handle three things at once: the K-element coupling matrix, +* the .tran "uic" flag, and a stiff multi-time-scale transient. +V1 inp 0 DC 1 +R10 inp x 0.5 +C2 x 0 2 +L8 x c 2H ic=2 +C12 x y 5 ic=2 +L9 y c 2H +C3 y 0 4 ic=5 +R11 y out 0.25 +R5 out 0 2 +I1 0 out DC 1 +L6 c n1 4H +R4 n1 0 1 +K1 L6 L8 -0.3535534 +K2 L6 L9 -0.3535534 +K3 L8 L9 -0.5 +.tran 0.2 200 0 0.2 uic`, }, { integrationMethod: 'gear2', diff --git a/packages/core/src/analysis/transient-driver.ts b/packages/core/src/analysis/transient-driver.ts index b0404a7..27cc68c 100644 --- a/packages/core/src/analysis/transient-driver.ts +++ b/packages/core/src/analysis/transient-driver.ts @@ -104,6 +104,13 @@ interface InternalTransientConfig { maxTimestep: number; /** Optional pre-computed DC solution. When provided, skips internal DC op point. */ initialSolution?: Float64Array; + /** + * The seed is a UIC initial-condition solution rather than a DC operating + * point. The companion model has no valid history at t=0, so the first + * step must bootstrap with Backward Euler and a cut dt — same pattern the + * driver uses after crossing a waveform breakpoint. + */ + seedIsUIC?: boolean; } class TransientSimImpl implements TransientSim { @@ -136,7 +143,15 @@ class TransientSimImpl implements TransientSim { if (config.initialSolution) { // Caller already computed DC — skip internal DC and seed directly. this.assembler.solution.set(config.initialSolution); - this.stampPrevB(); + if (config.seedIsUIC) { + // Force BE bootstrap on the first step: leave prevB undefined and + // trigger the post-breakpoint dt cut. The seeded UIC solution does + // not satisfy any prior companion equation, so any trap/gear-2 + // history term we built here would be mathematically inconsistent. + this.justCrossedBreakpoint = true; + } else { + this.stampPrevB(); + } } else { this.initDC(); } @@ -369,6 +384,7 @@ export function createDriverFromCompiled( timestep: number; maxTimestep: number; initialSolution?: Float64Array; + seedIsUIC?: boolean; }, ): TransientSim & { peekInitialStep(): TransientStep } { const impl = new TransientSimImpl(compiled, options, { @@ -376,6 +392,7 @@ export function createDriverFromCompiled( timestep: config.timestep, maxTimestep: config.maxTimestep, initialSolution: config.initialSolution, + seedIsUIC: config.seedIsUIC, }); return impl; } diff --git a/packages/core/src/analysis/transient.ts b/packages/core/src/analysis/transient.ts index b3099ba..c8f8705 100644 --- a/packages/core/src/analysis/transient.ts +++ b/packages/core/src/analysis/transient.ts @@ -12,6 +12,7 @@ export function solveTransient( analysis: TransientAnalysis, options: ResolvedOptions, initialSolution?: Float64Array, + seedIsUIC?: boolean, ): TransientResult { const { nodeNames, branchNames } = compiled; const driver = createDriverFromCompiled(compiled, options, { @@ -19,6 +20,7 @@ export function solveTransient( timestep: analysis.timestep, maxTimestep: analysis.maxTimestep ?? Math.min(analysis.timestep, analysis.stopTime / 50), initialSolution, + seedIsUIC, }); const timePoints: number[] = [0]; diff --git a/packages/core/src/analysis/uic.ts b/packages/core/src/analysis/uic.ts new file mode 100644 index 0000000..2029c67 --- /dev/null +++ b/packages/core/src/analysis/uic.ts @@ -0,0 +1,73 @@ +import type { Circuit, CompiledCircuit } from '../circuit.js'; +import type { ResolvedOptions } from '../types.js'; +import { solveDCOperatingPoint } from './dc.js'; +import { Inductor } from '../devices/inductor.js'; + +/** + * Compute the initial solution vector for a `.tran ... uic` run. + * + * Substitutes capacitors with `ic=` for DC voltage sources, inductors with + * `ic=` for DC current sources, and drops mutual-inductance K-elements + * (which contribute nothing at DC). Solves the resulting linear DC system, + * then maps node voltages and branch currents back to the *original* + * compiled circuit's solution layout — which uses different branch indices + * because inductors-with-ic become branchless current sources, and + * capacitors-with-ic gain new branches. + * + * Storage elements without `ic` keep their default DC behaviour: caps act + * as opens, inductors as shorts. This matches the ngspice convention where + * the unspecified element settles to the value implied by the surrounding + * topology rather than being forced to zero. + * + * @returns Initial solution sized for the original `compiled` system. + */ +export function computeUICInitialSolution( + circuit: Circuit, + compiled: CompiledCircuit, + options: ResolvedOptions, +): Float64Array { + const icCompiled = circuit.compileForUIC(); + const { assembler: icAsm } = solveDCOperatingPoint(icCompiled, options); + + const seed = new Float64Array(compiled.nodeCount + compiled.branchCount); + + // Map node voltages by name. The IC system has the same node set as the + // original (substitutions only swap device types, not topology). + for (const [name, origIdx] of compiled.nodeIndexMap) { + if (origIdx < 0) continue; + const icIdx = icCompiled.nodeIndexMap.get(name); + if (icIdx !== undefined && icIdx >= 0) { + seed[origIdx] = icAsm.solution[icIdx]; + } + } + + // Map branch currents. Original branches come from V sources, inductors + // (all of them), and controlled sources (E/H). After substitution: + // - V sources keep their branches (currents transfer directly). + // - Inductors WITH ic became I sources (no branch in IC system); use the + // declared ic as the seed current, with the same dot/sign convention + // as the inductor's own KCL stamp (branch current enters n+). + // - Inductors WITHOUT ic stayed as inductors (DC short, branch present). + // - E/H controlled sources keep their branches. + for (let i = 0; i < compiled.branchNames.length; i++) { + const name = compiled.branchNames[i]; + const origBranchAbs = compiled.nodeCount + i; + + // Look up the original device to know what kind of element this branch + // belongs to. + const dev = compiled.devices.find(d => d.name === name); + + if (dev instanceof Inductor && dev.ic !== undefined) { + seed[origBranchAbs] = dev.ic; + continue; + } + + // Otherwise, find the same name in the IC system's branch layout. + const icBranchIdx = icCompiled.branchNames.indexOf(name); + if (icBranchIdx >= 0) { + seed[origBranchAbs] = icAsm.solution[icCompiled.nodeCount + icBranchIdx]; + } + } + + return seed; +} diff --git a/packages/core/src/circuit.ts b/packages/core/src/circuit.ts index dfeebb4..0e1a70f 100644 --- a/packages/core/src/circuit.ts +++ b/packages/core/src/circuit.ts @@ -7,6 +7,7 @@ import { VoltageSource } from './devices/voltage-source.js'; import { CurrentSource } from './devices/current-source.js'; import { Capacitor } from './devices/capacitor.js'; import { Inductor } from './devices/inductor.js'; +import { MutualInductor } from './devices/mutual-inductor.js'; import { Diode } from './devices/diode.js'; import { BJT } from './devices/bjt.js'; import { MOSFET } from './devices/mosfet.js'; @@ -60,6 +61,11 @@ interface DeviceDescriptor { modelName?: string; params?: Record; controlSource?: string; + /** Initial condition (volts for caps, amps for inductors). */ + ic?: number; + /** For 'K' descriptors: names of the two inductors being coupled. */ + coupledA?: string; + coupledB?: string; } /** @@ -130,10 +136,10 @@ export class Circuit { * @param nodeNeg - Negative terminal node * @param capacitance - Capacitance value in farads */ - addCapacitor(name: string, nodePos: string, nodeNeg: string, capacitance: number): void { + addCapacitor(name: string, nodePos: string, nodeNeg: string, capacitance: number, ic?: number): void { this.nodeSet.add(nodePos); this.nodeSet.add(nodeNeg); - this.descriptors.push({ type: 'C', name, nodes: [nodePos, nodeNeg], value: capacitance }); + this.descriptors.push({ type: 'C', name, nodes: [nodePos, nodeNeg], value: capacitance, ic }); } /** @@ -144,10 +150,30 @@ export class Circuit { * @param nodeNeg - Negative terminal node * @param inductance - Inductance value in henries */ - addInductor(name: string, nodePos: string, nodeNeg: string, inductance: number): void { + addInductor(name: string, nodePos: string, nodeNeg: string, inductance: number, ic?: number): void { this.nodeSet.add(nodePos); this.nodeSet.add(nodeNeg); - this.descriptors.push({ type: 'L', name, nodes: [nodePos, nodeNeg], value: inductance }); + this.descriptors.push({ type: 'L', name, nodes: [nodePos, nodeNeg], value: inductance, ic }); + } + + /** + * Couple two inductors via a coupling coefficient k (K-element). + * + * The mutual inductance is M = k * sqrt(L_a * L_b), with |k| ≤ 1 for a + * physically realizable coupling. Both inductors must already exist via + * {@link addInductor}; their dot convention is implicit in their declared + * node order. + * + * @param name - K-element name (e.g., `'K1'`) + * @param indA - First inductor name (e.g., `'L1'`) + * @param indB - Second inductor name (e.g., `'L2'`) + * @param k - Coupling coefficient (positive or negative) + */ + addInductorCoupling(name: string, indA: string, indB: string, k: number): void { + this.descriptors.push({ + type: 'K', name, nodes: [], + coupledA: indA, coupledB: indB, value: k, + }); } /** @@ -361,7 +387,7 @@ export class Circuit { */ addAnalysis(type: 'op'): void; addAnalysis(type: 'dc', params: { source: string; start: number; stop: number; step: number }): void; - addAnalysis(type: 'tran', params: { timestep: number; stopTime: number; startTime?: number; maxTimestep?: number }): void; + addAnalysis(type: 'tran', params: { timestep: number; stopTime: number; startTime?: number; maxTimestep?: number; uic?: boolean }): void; addAnalysis(type: 'ac', params: { variation: 'dec' | 'oct' | 'lin'; points: number; startFreq: number; stopFreq: number }): void; addAnalysis(type: string, params?: Record): void { switch (type) { @@ -378,13 +404,14 @@ export class Circuit { }); break; case 'tran': { - const tranCmd: { type: 'tran'; timestep: number; stopTime: number; startTime?: number; maxTimestep?: number } = { + const tranCmd: { type: 'tran'; timestep: number; stopTime: number; startTime?: number; maxTimestep?: number; uic?: boolean } = { type: 'tran', timestep: params!.timestep as number, stopTime: params!.stopTime as number, }; if (params?.startTime !== undefined) tranCmd.startTime = params.startTime as number; if (params?.maxTimestep !== undefined) tranCmd.maxTimestep = params.maxTimestep as number; + if (params?.uic) tranCmd.uic = true; this._analyses.push(tranCmd); break; } @@ -450,9 +477,47 @@ export class Circuit { * @throws Error if a referenced subcircuit or control source is undefined * @throws {@link CycleError} if subcircuit instances form a circular dependency */ + /** + * Compile a UIC seed circuit: capacitors with `ic=` are replaced by DC + * voltage sources of that value, and inductors with `ic=` are replaced + * by DC current sources. Storage elements without `ic` keep their default + * DC behaviour (caps open, inductors short). K-elements are dropped (they + * have no DC contribution). Used to compute the t=0 seed for `.tran ... uic`. + */ + compileForUIC(): CompiledCircuit { + return this.compileWith(d => this.toUICDescriptor(d)); + } + + private toUICDescriptor(desc: DeviceDescriptor): DeviceDescriptor | null { + if (desc.type === 'C' && desc.ic !== undefined) { + return { + type: 'V', name: desc.name, nodes: desc.nodes, + waveform: { type: 'dc', value: desc.ic }, + }; + } + if (desc.type === 'L' && desc.ic !== undefined) { + // Inductor convention: positive branch current flows from n+ to n-, + // i.e., LEAVES n+. SPICE current source convention: positive value + // INJECTS into n+. Flip the node order so an I source of +ic models + // an inductor with branch current +ic. + return { + type: 'I', name: desc.name, nodes: [desc.nodes[1], desc.nodes[0]], + waveform: { type: 'dc', value: desc.ic }, + }; + } + if (desc.type === 'K') return null; // no DC effect + return desc; + } + compile(): CompiledCircuit { + return this.compileWith(d => d); + } + + private compileWith(transform: (d: DeviceDescriptor) => DeviceDescriptor | null): CompiledCircuit { // Pre-expand subcircuit instances into flat device descriptors - const expandedDescriptors = this.expandAllSubcircuits(); + const expandedDescriptors = this.expandAllSubcircuits() + .map(transform) + .filter((d): d is DeviceDescriptor => d !== null); // Collect all nodes from expanded descriptors for (const desc of expandedDescriptors) { @@ -503,12 +568,24 @@ export class Circuit { devices.push(new CurrentSource(desc.name, nodeIndices, resolveWaveform(desc.waveform))); break; case 'C': - devices.push(new Capacitor(desc.name, nodeIndices, desc.value!)); + devices.push(new Capacitor(desc.name, nodeIndices, desc.value!, desc.ic)); break; case 'L': { const bi = branchIndex++; branchNames.push(desc.name); - devices.push(new Inductor(desc.name, nodeIndices, bi, desc.value!)); + devices.push(new Inductor(desc.name, nodeIndices, bi, desc.value!, desc.ic)); + break; + } + case 'K': { + // Resolve the two coupled inductor names to their device instances. + const a = deviceMap.get(desc.coupledA!); + const b = deviceMap.get(desc.coupledB!); + if (!(a instanceof Inductor) || !(b instanceof Inductor)) { + throw new Error( + `K-element '${desc.name}' references unknown or non-inductor device(s): ${desc.coupledA}, ${desc.coupledB}`, + ); + } + devices.push(new MutualInductor(desc.name, a, b, desc.value!)); break; } case 'D': { diff --git a/packages/core/src/devices/capacitor.ts b/packages/core/src/devices/capacitor.ts index d6a7217..b72a885 100644 --- a/packages/core/src/devices/capacitor.ts +++ b/packages/core/src/devices/capacitor.ts @@ -8,6 +8,7 @@ export class Capacitor implements DeviceModel { readonly name: string, readonly nodes: number[], public capacitance: number, + public ic?: number, ) {} setParameter(value: number): void { diff --git a/packages/core/src/devices/index.ts b/packages/core/src/devices/index.ts index 45095a9..62f879a 100644 --- a/packages/core/src/devices/index.ts +++ b/packages/core/src/devices/index.ts @@ -2,6 +2,7 @@ export type { DeviceModel, StampContext } from './device.js'; export { Resistor } from './resistor.js'; export { Capacitor } from './capacitor.js'; export { Inductor } from './inductor.js'; +export { MutualInductor } from './mutual-inductor.js'; export { VoltageSource } from './voltage-source.js'; export { CurrentSource } from './current-source.js'; export { Diode } from './diode.js'; diff --git a/packages/core/src/devices/inductor.ts b/packages/core/src/devices/inductor.ts index 17576aa..81915fb 100644 --- a/packages/core/src/devices/inductor.ts +++ b/packages/core/src/devices/inductor.ts @@ -9,6 +9,7 @@ export class Inductor implements DeviceModel { readonly nodes: number[], readonly branchIndex: number, public inductance: number, + public ic?: number, ) { this.branches = [branchIndex]; } diff --git a/packages/core/src/devices/mutual-inductor.test.ts b/packages/core/src/devices/mutual-inductor.test.ts new file mode 100644 index 0000000..1197619 --- /dev/null +++ b/packages/core/src/devices/mutual-inductor.test.ts @@ -0,0 +1,75 @@ +import { describe, it, expect } from 'vitest'; +import { simulate } from '../simulate.js'; + +describe('Mutual inductance (K-element)', () => { + it('parses K-line and produces coupled-inductor response (Chua & Lin 8-7)', async () => { + // Full Chua & Lin page 343 benchmark — three coupled inductors with + // K1, K2, K3 and an asymmetric IC seed. Verifies the K-element + // contribution actually reaches the dynamic system: coupled vs. + // uncoupled traces differ by a measurable amount on the inductor + // currents. + const tplt = (kLines: string) => ` + V1 inp 0 DC 1 + R10 inp x 0.5 + C2 x 0 2 + L8 x c 2H ic=2 + C12 x y 5 ic=2 + L9 y c 2H + C3 y 0 4 ic=5 + R11 y out 0.25 + R5 out 0 2 + I1 0 out DC 1 + L6 c n1 4H + R4 n1 0 1 + ${kLines} + .tran 0.5 5 0 0.5 uic`; + const coupled = (await simulate(tplt( + 'K1 L6 L8 -0.3535534\nK2 L6 L9 -0.3535534\nK3 L8 L9 -0.5', + ))).transient!; + const uncoupled = (await simulate(tplt('* no coupling'))).transient!; + + const ic = coupled.current('L9'); + const iu = uncoupled.current('L9'); + let maxDiff = 0; + for (let i = 0; i < Math.min(ic.length, iu.length); i++) { + maxDiff = Math.max(maxDiff, Math.abs(ic[i] - iu[i])); + } + expect(maxDiff).toBeGreaterThan(0.1); + }); + + it('throws on K referencing an unknown inductor', async () => { + await expect(simulate(` + V1 a 0 DC 1 + L1 a 0 1m + K1 L1 L99 0.5 + .tran 1u 1m + `)).rejects.toThrow(/K-element/); + }); +}); + +describe('Initial conditions (.tran ... uic)', () => { + it('seeds capacitor voltages and inductor currents from ic= values', async () => { + // RC discharge: C with ic=5V, no source. With UIC, V across cap starts + // at 5V and decays through R toward 0. + const r = await simulate(` + C1 a 0 1u ic=5 + R1 a 0 1k + .tran 100u 5m uic + `); + const v = r.transient!.voltage('a'); + expect(v[0]).toBeCloseTo(5, 1); + expect(v[v.length - 1]).toBeLessThan(0.1); // 5 RC settled + }); + + it('seeds inductor branch current from ic=', async () => { + // RL freewheel: L with ic=2A discharging through R. + const r = await simulate(` + L1 a 0 1m ic=2 + R1 a 0 100 + .tran 1u 50u uic + `); + const i = r.transient!.current('L1'); + expect(i[0]).toBeCloseTo(2, 1); + expect(Math.abs(i[i.length - 1])).toBeLessThan(0.2); // ~5 L/R settled + }); +}); diff --git a/packages/core/src/devices/mutual-inductor.ts b/packages/core/src/devices/mutual-inductor.ts new file mode 100644 index 0000000..1193a72 --- /dev/null +++ b/packages/core/src/devices/mutual-inductor.ts @@ -0,0 +1,43 @@ +import type { DeviceModel, StampContext } from './device.js'; +import type { Inductor } from './inductor.js'; + +/** + * Mutual inductance (K-element) — couples two inductors with a coefficient k. + * + * The mutual inductance is M = k * sqrt(L_a * L_b), where |k| ≤ 1 for + * physically realizable couplings. The branch equations of the two + * coupled inductors become: + * V_a = L_a * dI_a/dt + M * dI_b/dt + * V_b = M * dI_a/dt + L_b * dI_b/dt + * + * In MNA terms this adds two off-diagonal entries to the C matrix at + * (branch_a, branch_b) and (branch_b, branch_a), each equal to -M + * (matching the existing -L sign convention used by {@link Inductor}). + * + * The device contributes nothing at DC and nothing to the conductance + * matrix — only to the dynamic (C) matrix used by the companion model. + */ +export class MutualInductor implements DeviceModel { + readonly nodes: number[] = []; + readonly branches: number[] = []; + readonly isNonlinear = false; + + constructor( + readonly name: string, + readonly indA: Inductor, + readonly indB: Inductor, + public k: number, + ) {} + + stamp(_ctx: StampContext): void { + // No DC contribution; mutual coupling is purely a dynamic (dI/dt) effect. + } + + stampDynamic(ctx: StampContext): void { + const M = this.k * Math.sqrt(this.indA.inductance * this.indB.inductance); + const ba = ctx.numNodes + this.indA.branchIndex; + const bb = ctx.numNodes + this.indB.branchIndex; + ctx.stampC(ba, bb, -M); + ctx.stampC(bb, ba, -M); + } +} diff --git a/packages/core/src/parser/index.ts b/packages/core/src/parser/index.ts index e3c78d0..29d6d3b 100644 --- a/packages/core/src/parser/index.ts +++ b/packages/core/src/parser/index.ts @@ -3,6 +3,19 @@ import { ParseError } from '../errors.js'; import { tokenizeNetlist, parseNumber } from './tokenizer.js'; import { parseModelCard } from './model-parser.js'; import { parseSourceWaveform, parseInstanceParams } from './waveform-parser.js'; + +/** Extract `IC=N` from a token list (case-insensitive). Returns undefined if absent. */ +function pickIC(tokens: string[], startIdx: number): number | undefined { + for (let i = startIdx; i < tokens.length; i++) { + const t = tokens[i]; + const eq = t.indexOf('='); + if (eq <= 0) continue; + if (t.slice(0, eq).toLowerCase() === 'ic') { + return parseNumber(t.slice(eq + 1)); + } + } + return undefined; +} import { preprocess } from './preprocessor.js'; import type { IncludeResolver } from '../types.js'; @@ -137,11 +150,18 @@ function parseDotCommand(circuit: Circuit, tokens: string[], lineNumber: number) break; } case '.TRAN': { - const timestep = parseNumber(tokens[1]); - const stopTime = parseNumber(tokens[2]); - const startTime = tokens[3] ? parseNumber(tokens[3]) : undefined; - const maxTimestep = tokens[4] ? parseNumber(tokens[4]) : undefined; - circuit.addAnalysis('tran', { timestep, stopTime, startTime, maxTimestep }); + // Trailing 'UIC' flag is a positional keyword, not a number. + const numericTokens: string[] = []; + let uic = false; + for (let i = 1; i < tokens.length; i++) { + if (tokens[i].toUpperCase() === 'UIC') uic = true; + else numericTokens.push(tokens[i]); + } + const timestep = parseNumber(numericTokens[0]); + const stopTime = parseNumber(numericTokens[1]); + const startTime = numericTokens[2] ? parseNumber(numericTokens[2]) : undefined; + const maxTimestep = numericTokens[3] ? parseNumber(numericTokens[3]) : undefined; + circuit.addAnalysis('tran', { timestep, stopTime, startTime, maxTimestep, uic }); break; } case '.AC': { @@ -213,12 +233,20 @@ function parseDevice(circuit: Circuit, tokens: string[], lineNumber: number): vo } case 'C': { const value = parseNumber(tokens[3]); - circuit.addCapacitor(name, tokens[1], tokens[2], value); + const ic = pickIC(tokens, 4); + circuit.addCapacitor(name, tokens[1], tokens[2], value, ic); break; } case 'L': { const value = parseNumber(tokens[3]); - circuit.addInductor(name, tokens[1], tokens[2], value); + const ic = pickIC(tokens, 4); + circuit.addInductor(name, tokens[1], tokens[2], value, ic); + break; + } + case 'K': { + // K + const k = parseNumber(tokens[3]); + circuit.addInductorCoupling(name, tokens[1], tokens[2], k); break; } case 'V': { diff --git a/packages/core/src/parser/tokenizer.ts b/packages/core/src/parser/tokenizer.ts index 0fe4001..5891885 100644 --- a/packages/core/src/parser/tokenizer.ts +++ b/packages/core/src/parser/tokenizer.ts @@ -36,13 +36,25 @@ export function parseNumber(token: string): number { } // Single-char suffixes — case-sensitive (m = milli, M = mega) - const suffixMatch = trimmed.match(/^([+-]?[\d.]+(?:[eE][+-]?\d+)?)([TtGgKkMmUuNnPpFf])$/); + // Trailing alphabetic chars after the multiplier are treated as unit + // annotations (e.g. "1uF", "10kHz") and ignored — SPICE convention. + const suffixMatch = trimmed.match(/^([+-]?[\d.]+(?:[eE][+-]?\d+)?)([TtGgKkMmUuNnPpFf])([A-Za-z]*)$/); if (suffixMatch) { const exp = SI_SUFFIX_MAP[suffixMatch[2]]; const val = Number(suffixMatch[1] + exp); if (!isNaN(val)) return val; } + // Plain number followed by unit letters (e.g. "2H", "5V", "100Ohm"). + // Only matches when the entire trailing portion is alphabetic — keeps + // ambiguous cases like "5F" out of this branch since the SI suffix + // parse above will already have handled them as femto-units. + const plainWithUnit = trimmed.match(/^([+-]?(?:\d+\.?\d*|\.\d+)(?:[eE][+-]?\d+)?)([A-Za-z]+)$/); + if (plainWithUnit) { + const val = Number(plainWithUnit[1]); + if (!isNaN(val)) return val; + } + throw new Error(`Cannot parse number: '${token}'`); } diff --git a/packages/core/src/simulate.ts b/packages/core/src/simulate.ts index 5ddbfa1..ad39f2c 100644 --- a/packages/core/src/simulate.ts +++ b/packages/core/src/simulate.ts @@ -6,6 +6,7 @@ import type { TransientAnalysis, ACAnalysis, ResolvedOptions } from './types.js' import { resolveOptions } from './types.js'; import { solveDCOperatingPoint } from './analysis/dc.js'; import { solveTransient } from './analysis/transient.js'; +import { computeUICInitialSolution } from './analysis/uic.js'; import { solveAC } from './analysis/ac.js'; import { solveDCSweep } from './analysis/dc-sweep.js'; import { solveStep, generateStepValues } from './analysis/step.js'; @@ -87,8 +88,10 @@ export async function simulate( } case 'tran': { const opts = resolveOptions(options, analysis.stopTime); - const { assembler: dcAsm } = solveDCOperatingPoint(compiled, opts); - result.transient = solveTransient(compiled, analysis, opts, dcAsm.solution); + const seed = analysis.uic + ? computeUICInitialSolution(circuit, compiled, opts) + : solveDCOperatingPoint(compiled, opts).assembler.solution; + result.transient = solveTransient(compiled, analysis, opts, seed, analysis.uic); break; } case 'ac': { diff --git a/packages/core/src/types.ts b/packages/core/src/types.ts index c525d98..e462655 100644 --- a/packages/core/src/types.ts +++ b/packages/core/src/types.ts @@ -36,6 +36,12 @@ export interface TransientAnalysis { startTime?: number; /** Maximum allowed timestep in seconds */ maxTimestep?: number; + /** + * Use Initial Conditions: skip the DC operating point and seed the + * transient from `ic=` values declared on capacitors and inductors. + * Storage elements without an explicit `ic` settle from V=0 / I=0. + */ + uic?: boolean; } /** AC small-signal analysis (`.ac`). Frequency sweep. */