init

4 files changed, 263 insertions, 0 deletions

diff --git a/src/physics/electricField.test.ts b/src/physics/electricField.test.ts
new file mode 100644
index 0000000..2e01309
--- /dev/null
+++ b/src/physics/electricField.test.ts
@@ -0,0 +1,21 @@
+import { describe, expect, it } from 'vitest';
+import { electricField, COULOMB_CONSTANT } from './electricField';
+
+describe('Coulomb electric field', () => {
+  it('points away from a positive charge with inverse-square magnitude', () => {
+    const field = electricField({ x: 2, y: 0, z: 0 }, [{ x: 0, y: 0, z: 0, coulombs: 1e-9 }]);
+    expect(field.x).toBeCloseTo(COULOMB_CONSTANT * 1e-9 / 4);
+    expect(field.y).toBe(0);
+    expect(field.z).toBe(0);
+  });
+  it('superposes opposite charges in the same direction at the midpoint', () => {
+    const field = electricField({ x: 0, y: 0, z: 0 }, [
+      { x: -1, y: 0, z: 0, coulombs: 1e-9 },
+      { x: 1, y: 0, z: 0, coulombs: -1e-9 },
+    ]);
+    expect(field.x).toBeCloseTo(2 * COULOMB_CONSTANT * 1e-9);
+  });
+  it('omits singular samples within the visualization exclusion radius', () => {
+    expect(electricField({ x: 0, y: 0, z: 0 }, [{ x: 0, y: 0, z: 0, coulombs: 1e-9 }])).toEqual({ x: 0, y: 0, z: 0 });
+  });
+});
diff --git a/src/physics/electricField.ts b/src/physics/electricField.ts
new file mode 100644
index 0000000..26975e3
--- /dev/null
+++ b/src/physics/electricField.ts
@@ -0,0 +1,22 @@
+export interface Point { x: number; y: number; z: number }
+export interface Charge extends Point { coulombs: number }
+export const COULOMB_CONSTANT = 8.9875517923e9;
+
+/** SI units: positions in metres, charges in coulombs, output in N/C.
+ * Samples inside 0.18 m of a charge are omitted to avoid visual singularities.
+ */
+export function electricField(point: Point, charges: readonly Charge[]): Point {
+  const field = { x: 0, y: 0, z: 0 };
+  for (const charge of charges) {
+    const dx = point.x - charge.x;
+    const dy = point.y - charge.y;
+    const dz = point.z - charge.z;
+    const r2 = dx * dx + dy * dy + dz * dz;
+    if (r2 < 0.18 ** 2) continue;
+    const factor = COULOMB_CONSTANT * charge.coulombs / (r2 * Math.sqrt(r2));
+    field.x += factor * dx;
+    field.y += factor * dy;
+    field.z += factor * dz;
+  }
+  return field;
+}
diff --git a/src/physics/models.test.ts b/src/physics/models.test.ts
new file mode 100644
index 0000000..3ff6487
--- /dev/null
+++ b/src/physics/models.test.ts
@@ -0,0 +1,108 @@
+import { describe, expect, it } from 'vitest';
+import { electricField } from './electricField';
+import { potential, forceOn, pairEnergy, gaussianFlux, capacitor, conductor, semiconductor } from './models';
+
+// Fixtures use hand-calculated SI results, independent of model implementations.
+describe('charge, potential and energy', () => {
+  const source = [{ x: 0, y: 0, z: 0, coulombs: 2e-9 }];
+  it('gives scalar potential tending to zero at infinity', () => {
+    expect(potential({ x: 2, y: 0, z: 0 }, source)).toBeCloseTo(8.98755, 4);
+    expect(potential({ x: 1e9, y: 0, z: 0 }, source)).toBeCloseTo(0, 6);
+  });
+  it('makes the field the negative spatial gradient of potential', () => {
+    const h = 1e-4;
+    const numerical = -(potential({ x: 2 + h, y: 0, z: 0 }, source) - potential({ x: 2 - h, y: 0, z: 0 }, source)) / (2 * h);
+    expect(electricField({ x: 2, y: 0, z: 0 }, source).x).toBeCloseTo(numerical, 5);
+  });
+  it('reverses force for a negative test charge', () => {
+    expect(forceOn({ x: 2, y: 0, z: 0, coulombs: -3e-9 }, source).x).toBeCloseTo(-13.48132768845e-9, 16);
+  });
+  it('gives negative interaction energy for opposite charges', () => {
+    expect(pairEnergy(2e-9, -3e-9, 2)).toBeCloseTo(-26.9626553769e-9, 16);
+  });
+  it('marks potential at a charge as undefined instead of inventing a finite value', () => {
+    expect(potential({ x: 0, y: 0, z: 0 }, source)).toBeNaN();
+    expect(() => pairEnergy(1, 1, 0)).toThrow(RangeError);
+  });
+});
+
+describe('Gauss law', () => {
+  it('counts only enclosed charge, including its sign', () => {
+    const charges = [{ x: 0, y: 0, z: 0, coulombs: -2e-9 }, { x: 3, y: 0, z: 0, coulombs: 5e-9 }];
+    expect(gaussianFlux(charges, 2)).toBeCloseTo(-225.8818, 3);
+    expect(gaussianFlux(charges, 4)).toBeCloseTo(338.8227, 3);
+  });
+  it('does not depend on surface size when enclosed charge is unchanged', () => {
+    const charges = [{ x: 0, y: 0, z: 0, coulombs: 1e-9 }];
+    expect(gaussianFlux(charges, 1)).toBe(gaussianFlux(charges, 2));
+  });
+  it('marks charge on the Gaussian surface as undefined', () => {
+    expect(gaussianFlux([{ x: 2, y: 0, z: 0, coulombs: 1e-9 }], 2)).toBeNaN();
+  });
+});
+
+describe('parallel plate capacitor', () => {
+  it('converts physical geometry into C, Q, E and stored energy', () => {
+    const c = capacitor({ area: 0.01, gap: 0.001, relativePermittivity: 1, voltage: 10 });
+    expect(c.capacitance).toBeCloseTo(8.8541878128e-11, 18);
+    expect(c.charge).toBeCloseTo(8.8541878128e-10, 18);
+    expect(c.field).toBe(10000);
+    expect(c.energy).toBeCloseTo(4.4270939064e-9, 17);
+  });
+  it('increases free charge at fixed voltage when dielectric is inserted', () => {
+    const c = capacitor({ area: 0.01, gap: 0.001, relativePermittivity: 4, voltage: 10 });
+    expect(c.field).toBe(10000);
+    expect(c.charge).toBeCloseTo(3.54167512512e-9, 17);
+    expect(c.polarization).toBeCloseTo(2.65625634384e-7, 15);
+  });
+  it('reduces voltage and energy at fixed free charge', () => {
+    const c = capacitor({ area: 0.01, gap: 0.001, relativePermittivity: 4, voltage: 10, isolated: true });
+    expect(c.charge).toBeCloseTo(8.8541878128e-10, 18);
+    expect(c.voltage).toBeCloseTo(2.5);
+    expect(c.field).toBeCloseTo(2500);
+    expect(c.energy).toBeCloseTo(1.1067734766e-9, 17);
+  });
+  it('rejects impossible geometry', () => {
+    expect(() => capacitor({ area: 0.01, gap: 0, relativePermittivity: 1, voltage: 10 })).toThrow(RangeError);
+  });
+});
+
+describe('current and drift', () => {
+  it('follows Ohm law and sends electrons opposite conventional current', () => {
+    const wire = conductor({ voltage: 2, resistivity: 1.68e-8, length: 1, area: 1e-6, carrierDensity: 8.5e28 });
+    expect(wire.resistance).toBeCloseTo(0.0168, 6);
+    expect(wire.current).toBeCloseTo(119.047619, 5);
+    expect(wire.electronDrift).toBeCloseTo(-0.00874160935, 10);
+    expect(wire.power).toBeCloseTo(238.095238, 5);
+  });
+  it('reverses current but keeps dissipated power positive with reversed voltage', () => {
+    const wire = conductor({ voltage: -2, resistivity: 1.68e-8, length: 1, area: 1e-6, carrierDensity: 8.5e28 });
+    expect(wire.current).toBeLessThan(0);
+    expect(wire.electronDrift).toBeGreaterThan(0);
+    expect(wire.power).toBeGreaterThan(0);
+  });
+});
+
+describe('semiconductor equilibrium at 300 K', () => {
+  it('has equal intrinsic electrons and holes', () => {
+    const s = semiconductor({ intrinsicDensity: 1e16, netDoping: 0, electronMobility: 0.135, holeMobility: 0.048 });
+    expect(s.electrons).toBe(1e16);
+    expect(s.holes).toBe(1e16);
+    expect(s.fermiShift).toBe(0);
+    expect(s.conductivity).toBeCloseTo(0.0002931983, 9);
+  });
+  it('conserves mass action and charge neutrality even at heavy donor doping', () => {
+    const s = semiconductor({ intrinsicDensity: 1e16, netDoping: 1e22, electronMobility: 0.135, holeMobility: 0.048 });
+    expect(s.electrons / 1e22).toBeCloseTo(1, 10);
+    expect(s.holes).toBeCloseTo(1e10, 0);
+    expect(s.electrons * s.holes / 1e32).toBeCloseTo(1, 12);
+    expect((s.electrons - s.holes) / 1e22).toBeCloseTo(1, 12);
+    expect(s.fermiShift).toBeCloseTo(0.3571586, 5);
+  });
+  it('makes holes majority carriers for acceptor doping', () => {
+    const s = semiconductor({ intrinsicDensity: 1e16, netDoping: -1e22, electronMobility: 0.135, holeMobility: 0.048 });
+    expect(s.holes / 1e22).toBeCloseTo(1, 10);
+    expect(s.electrons).toBeCloseTo(1e10, 0);
+    expect(s.fermiShift).toBeLessThan(0);
+  });
+});
diff --git a/src/physics/models.ts b/src/physics/models.ts
new file mode 100644
index 0000000..8fb121f
--- /dev/null
+++ b/src/physics/models.ts
@@ -0,0 +1,112 @@
+import { COULOMB_CONSTANT, electricField, type Charge, type Point } from './electricField';
+
+export const VACUUM_PERMITTIVITY = 8.8541878128e-12;
+export const ELEMENTARY_CHARGE = 1.602176634e-19;
+const THERMAL_ENERGY_EV = 8.617333262e-5 * 300;
+
+function positive(...values: number[]) {
+  if (values.some(value => !Number.isFinite(value) || value <= 0)) {
+    throw new RangeError('Physical dimensions and material constants must be finite and positive.');
+  }
+}
+
+/** Vacuum potential in volts, referenced to infinity. Undefined inside the
+ * same 0.18 m visualization exclusion used by electricField. */
+export function potential(point: Point, charges: readonly Charge[]): number {
+  let volts = 0;
+  for (const charge of charges) {
+    if (charge.coulombs === 0) continue;
+    const r = Math.hypot(point.x - charge.x, point.y - charge.y, point.z - charge.z);
+    if (r < 0.18) return NaN;
+    volts += COULOMB_CONSTANT * charge.coulombs / r;
+  }
+  return volts;
+}
+
+/** Force on a test charge, excluding that charge from the supplied sources. */
+export function forceOn(test: Charge, sources: readonly Charge[]): Point {
+  const field = electricField(test, sources);
+  return { x: field.x * test.coulombs, y: field.y * test.coulombs, z: field.z * test.coulombs };
+}
+
+export function pairEnergy(q1: number, q2: number, separation: number): number {
+  positive(separation);
+  return COULOMB_CONSTANT * q1 * q2 / separation;
+}
+
+/** Exact signed flux for an origin-centred spherical Gaussian surface.
+ * A point charge exactly on the surface makes the integral undefined. */
+export function gaussianFlux(charges: readonly Charge[], radius: number): number {
+  positive(radius);
+  let enclosed = 0;
+  for (const charge of charges) {
+    if (charge.coulombs === 0) continue;
+    const r = Math.hypot(charge.x, charge.y, charge.z);
+    if (Math.abs(r - radius) < 1e-9) return NaN;
+    if (r < radius) enclosed += charge.coulombs;
+  }
+  return enclosed / VACUUM_PERMITTIVITY;
+}
+
+export interface CapacitorParameters {
+  area: number;
+  gap: number;
+  relativePermittivity: number;
+  voltage: number;
+  /** If isolated, voltage means the initial vacuum charging voltage. */
+  isolated?: boolean;
+}
+export function capacitor(p: CapacitorParameters) {
+  positive(p.area, p.gap, p.relativePermittivity);
+  const vacuumCapacitance = VACUUM_PERMITTIVITY * p.area / p.gap;
+  const capacitance = vacuumCapacitance * p.relativePermittivity;
+  const voltage = p.isolated ? p.voltage / p.relativePermittivity : p.voltage;
+  const charge = capacitance * voltage;
+  const field = voltage / p.gap;
+  return {
+    capacitance, voltage, charge, field,
+    energy: 0.5 * capacitance * voltage ** 2,
+    polarization: VACUUM_PERMITTIVITY * (p.relativePermittivity - 1) * field,
+  };
+}
+
+export interface ConductorParameters {
+  voltage: number;
+  resistivity: number;
+  length: number;
+  area: number;
+  carrierDensity: number;
+}
+export function conductor(p: ConductorParameters) {
+  positive(p.resistivity, p.length, p.area, p.carrierDensity);
+  const resistance = p.resistivity * p.length / p.area;
+  const current = p.voltage / resistance;
+  return {
+    resistance, current,
+    power: p.voltage * current,
+    density: current / p.area,
+    electronDrift: -current / (p.carrierDensity * ELEMENTARY_CHARGE * p.area),
+  };
+}
+
+export interface SemiconductorParameters {
+  intrinsicDensity: number;
+  /** Fully ionized donors minus acceptors, in m^-3. */
+  netDoping: number;
+  electronMobility: number;
+  holeMobility: number;
+}
+/** Non-degenerate thermal equilibrium at 300 K. Solving for the minority
+ * through mass action avoids catastrophic subtraction at large doping. */
+export function semiconductor(p: SemiconductorParameters) {
+  positive(p.intrinsicDensity, p.electronMobility, p.holeMobility);
+  const majority = (Math.abs(p.netDoping) + Math.hypot(p.netDoping, 2 * p.intrinsicDensity)) / 2;
+  const minority = p.intrinsicDensity ** 2 / majority;
+  const electrons = p.netDoping >= 0 ? majority : minority;
+  const holes = p.netDoping >= 0 ? minority : majority;
+  return {
+    electrons, holes,
+    fermiShift: THERMAL_ENERGY_EV * Math.log(electrons / p.intrinsicDensity),
+    conductivity: ELEMENTARY_CHARGE * (electrons * p.electronMobility + holes * p.holeMobility),
+  };
+}