luke
/
ngspice
mirror of https://git.code.sf.net/p/ngspice/ngspice
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

pa-140

This commit is contained in:
Meisam 2026-07-28 17:23:48 +02:00 committed by Holger Vogt
parent d3d4237313
commit 8cc60336ec
8 changed files with 467 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
src/frontend/Makefile.am
View File

@ -59,6 +59,8 @@ libfte_la_SOURCES = \
com_qpac.h \
com_qpnoise.c \
com_qpnoise.h \
com_hbosc.c \
com_hbosc.h \
com_hb.c \
com_hb.h \
com_checkpoint.c \

2
src/frontend/com_commands.h
View File

@ -10,6 +10,8 @@ void com_optimize(wordlist *wl); /* Enhancement-130 */
void com_qpss(wordlist *wl); /* Enhancement-133 */
void com_qpac(wordlist *wl); /* Enhancement-137 */
void com_qpnoise(wordlist *wl); /* Enhancement-138 */
void com_hbosc(wordlist *wl); /* Enhancement-140 */
void com_phasenoise(wordlist *wl); /* Enhancement-140 */
void com_hb(wordlist *wl); /* Enhancement-134 */
void com_savestate(wordlist *wl); /* Enhancement-131 */
void com_loadstate(wordlist *wl); /* Enhancement-131 */

173
src/frontend/com_hbosc.c Normal file
View File

@ -0,0 +1,173 @@
/**********
Enhancement-140: autonomous harmonic balance for oscillators + phase noise.
hbosc <oscnode> <K> [fguess] [tstab] -- autonomous HB: find the oscillator's
steady state (harmonics + frequency)
phasenoise <fstart> <fstop> [points] -- the phase-noise spectrum L(df)
`hbosc` runs a short transient (from the deck's .ic) to seed the limit cycle, estimates
the oscillation frequency and amplitude, and hands them to HBOSCanalyze() which refines
(V, w0) by a bordered Newton and retains the operating point. `phasenoise` then extracts
the perturbation projection vector (PPV) and folds the device noise to L(df). The engines
are in spicelib/analysis/dcpss.c.
**********/
#include "ngspice/ngspice.h"
#include "ngspice/cpdefs.h"
#include "ngspice/cktdefs.h"
#include "ngspice/ftedefs.h"
#include "ngspice/dvec.h"
#include "ngspice/fteext.h"
#include "ngspice/wordlist.h"
#include "ngspice/cpextern.h"
#include "circuits.h"
#include "com_hbosc.h"
static void hbosc_run_cmd(const char *cmdstr)
{
wordlist *wl = cp_lexer((char *) cmdstr);
int i;
if (!wl || !wl->wl_word) { if (wl) wl_free(wl); return; }
for (i = 0; cp_coms[i].co_comname; i++)
if (strcasecmp(cp_coms[i].co_comname, wl->wl_word) == 0)
break;
if (cp_coms[i].co_comname && cp_coms[i].co_func)
cp_coms[i].co_func(wl->wl_next);
wl_free(wl);
}
static double hboscnum(const char *w)
{
char *s = (char *) w;
double v = 0.0;
if (ft_numparse(&s, FALSE, &v) < 0)
v = atof(w);
return v;
}
static int hbosc_node(CKTcircuit *ckt, const char *name)
{
int numNames = 0, i, num = 0;
IFuid *nameList = NULL;
if (CKTnames(ckt, &numNames, &nameList) != OK || !nameList)
return 0;
for (i = 0; i < numNames; i++)
if (nameList[i] && strcmp((const char *) nameList[i], name) == 0) { num = i + 1; break; }
tfree(nameList);
return num;
}
void
com_hbosc(wordlist *wl)
{
CKTcircuit *ckt;
const char *oscname;
double fguess = 0.0, tstab = 0.0, ampseed, f0est;
int oscNode, K, verbose, err;
char cmd[128];
struct pnode *pn;
struct dvec *v, *sc;
double *tt, *vv;
int n, i0, i, ncross;
double vmax, tlast, tprev;
if (!ft_curckt || !ft_curckt->ci_ckt) {
fprintf(cp_err, "Error: hbosc: there is no circuit loaded.\n");
return;
}
ckt = ft_curckt->ci_ckt;
if (!wl || !wl->wl_next) {
fprintf(cp_err, "Usage: hbosc <oscnode> <K> [fguess] [tstab] (the deck needs a .ic to start the oscillation)\n");
return;
}
oscname = wl->wl_word;
K = (int) hboscnum(wl->wl_next->wl_word);
if (wl->wl_next->wl_next) {
fguess = hboscnum(wl->wl_next->wl_next->wl_word);
if (wl->wl_next->wl_next->wl_next)
tstab = hboscnum(wl->wl_next->wl_next->wl_next->wl_word);
}
if (K < 1) { fprintf(cp_err, "Error: hbosc: need K >= 1.\n"); return; }
if (fguess <= 0.0) fguess = 1e6; /* a default if none given */
if (tstab <= 0.0) tstab = 300.0 / fguess; /* ~300 periods to settle */
/* build up the limit cycle with a transient (uic honours the deck's .ic) */
(void) snprintf(cmd, sizeof cmd, "tran %.6g %.6g uic", 1.0/(fguess*40.0), tstab);
hbosc_run_cmd(cmd);
/* fetch v(oscnode) + its time scale */
(void) snprintf(cmd, sizeof cmd, "v(%s)", oscname);
pn = ft_getpnames_from_string(cmd, TRUE);
if (!pn) { fprintf(cp_err, "Error: hbosc: cannot read v(%s).\n", oscname); return; }
v = ft_evaluate(pn);
sc = (v && v->v_scale) ? v->v_scale : vec_get("time");
if (!v || !isreal(v) || v->v_length < 8 || !sc || sc->v_length < v->v_length) {
fprintf(cp_err, "Error: hbosc: no usable transient for v(%s).\n", oscname);
if (pn && !pn->pn_value && v) vec_free(v);
if (pn) free_pnode(pn);
return;
}
tt = sc->v_realdata; vv = v->v_realdata; n = v->v_length;
/* estimate the oscillation: amplitude = max|v| over the last third; frequency from
* the mean spacing of upward zero crossings there. */
i0 = (2*n)/3;
vmax = 0.0;
for (i = i0; i < n; i++) if (fabs(vv[i]) > vmax) vmax = fabs(vv[i]);
ncross = 0; tlast = tprev = 0.0;
for (i = i0 + 1; i < n; i++)
if (vv[i-1] <= 0.0 && vv[i] > 0.0) { /* upward crossing */
double tc = tt[i-1] + (tt[i]-tt[i-1]) * (-vv[i-1])/(vv[i]-vv[i-1]);
if (ncross == 0) tprev = tc;
tlast = tc; ncross++;
}
f0est = (ncross > 1 && tlast > tprev) ? (double)(ncross-1)/(tlast-tprev) : fguess;
ampseed = vmax;
oscNode = hbosc_node(ckt, oscname);
if (pn && !pn->pn_value && v) vec_free(v);
if (pn) free_pnode(pn);
if (oscNode <= 0) { fprintf(cp_err, "Error: hbosc: unknown oscillator node '%s'.\n", oscname); return; }
if (vmax < 1e-9) {
fprintf(cp_err, "Error: hbosc: no oscillation detected (add a `.ic` to start it).\n");
return;
}
verbose = cp_getvar("hbosc_verbose", CP_BOOL, NULL, 0);
err = HBOSCanalyze(ckt, oscNode, K, 0, f0est, ampseed, 60, 1e-11, verbose ? 1 : 0);
if (err != OK)
fprintf(cp_err, "hbosc: autonomous harmonic balance did not complete (error %d).\n", err);
}
void
com_phasenoise(wordlist *wl)
{
CKTcircuit *ckt;
double fstart, fstop;
int npts = 21, verbose, err;
if (!ft_curckt || !ft_curckt->ci_ckt) {
fprintf(cp_err, "Error: phasenoise: there is no circuit loaded.\n");
return;
}
ckt = ft_curckt->ci_ckt;
if (!wl || !wl->wl_next) {
fprintf(cp_err, "Usage: phasenoise <fstart> <fstop> [points] (run `hbosc` first)\n");
return;
}
fstart = hboscnum(wl->wl_word);
fstop = hboscnum(wl->wl_next->wl_word);
if (wl->wl_next->wl_next)
npts = (int) hboscnum(wl->wl_next->wl_next->wl_word);
if (fstart <= 0.0 || fstop < fstart) {
fprintf(cp_err, "Error: phasenoise: need 0 < fstart <= fstop.\n");
return;
}
if (npts < 1) npts = 1;
verbose = cp_getvar("phasenoise_verbose", CP_BOOL, NULL, 0);
err = PhaseNoiseAnalyze(ckt, fstart, fstop, npts, verbose ? 1 : 0);
if (err != OK)
fprintf(cp_err, "phasenoise: did not complete (error %d).\n", err);
}

8
src/frontend/com_hbosc.h Normal file
View File

@ -0,0 +1,8 @@
#ifndef ngspice_COM_HBOSC_H
#define ngspice_COM_HBOSC_H
/* Enhancement-140: autonomous HB for oscillators + phase noise. */
void com_hbosc(wordlist *wl);
void com_phasenoise(wordlist *wl);
#endif

8
src/frontend/commands.c
View File

@ -444,6 +444,14 @@ struct comm spcp_coms[] = {
{ 040, 040, 040, 040 }, E_DEFHMASK, 2, LOTS,
NULL,
"output_node f_in [cyclo] : two-tone QPnoise -- output/input noise density at f_in, folding device noise over all sidebands around the `qpss ... hb` operating point (`cyclo` = cyclostationary PSD)." },
{ "hbosc", com_hbosc, TRUE, FALSE, /* Enhancement-140 */
{ 040, 040, 040, 040 }, E_DEFHMASK, 2, LOTS,
NULL,
"oscnode K [fguess] [tstab] : autonomous harmonic balance -- an oscillator's steady state (harmonics + oscillation frequency), seeded from a transient." },
{ "phasenoise", com_phasenoise, TRUE, FALSE, /* Enhancement-140 */
{ 040, 040, 040, 040 }, E_DEFHMASK, 2, LOTS,
NULL,
"fstart fstop [points] : oscillator phase-noise spectrum L(df) via the PPV, around the `hbosc` operating point." },
{ "hb", com_hb, TRUE, FALSE, /* Enhancement-134 */
{ 040, 040, 040, 040 }, E_DEFHMASK, 2, LOTS,
NULL,

2
src/include/ngspice/cktdefs.h
View File

@ -479,6 +479,8 @@ extern int HBanalyze(CKTcircuit *, double f0, int K, int P, int maxiter, double
extern int QPSShb(CKTcircuit *, double f1, double f2, int K1, int K2, int P1, int P2, int maxiter, double tol, int verbose); /* E-136 */
extern int QPACanalyze(CKTcircuit *, double f_in, int verbose); /* E-137 */
extern int QPnoiseAnalyze(CKTcircuit *, int outNode, double f_in, int cyclo, int verbose); /* E-138 / -139 */
extern int HBOSCanalyze(CKTcircuit *, int oscNode, int K, int P, double f0seed, double ampseed, int maxiter, double tol, int verbose); /* E-140 */
extern int PhaseNoiseAnalyze(CKTcircuit *, double fstart, double fstop, int npts, int verbose); /* E-140 */
#endif

270
src/spicelib/analysis/dcpss.c
View File

@ -2425,6 +2425,276 @@ QPnoiseAnalyze(CKTcircuit *ckt, int outNode, double f_in, int cyclo, int verbose
}
/* ======================================================================
* Enhancement-140: oscillator phase noise.
*
* (1) HBOSCanalyze -- AUTONOMOUS harmonic balance. An oscillator has no
* driving source, so the HB residual F(V) = I_R + [dq/dt] = 0 is solved
* for the harmonics V *and* the unknown oscillation frequency w0, with a
* phase gauge (the solution's phase is free). The Jacobian dF/dV is the
* conversion matrix H -- SINGULAR, its right null space the phase mode
* u_k = jk V_k (shifting the oscillator's phase is a symmetry). Newton on
* the bordered system [ H dF/dw0 ; u*^T 0 ] [dV; dw0] = [-F; 0]
* (nonsingular by bordering) refines (V, w0) from a transient seed.
* (2) PhaseNoiseAnalyze -- the phase-noise spectrum L(df). The PPV (Demir's
* perturbation projection vector) is the LEFT null vector of H, normalized
* so <v1, phase-mode> = 1. The phase diffusion constant c is the device
* noise folded through the PPV (same machinery as pnoise, transfer -> PPV),
* and L(df) = 10 log10( f0^2 c / df^2 ) -- the classic 1/df^2 (-20 dB/dec)
* oscillator phase noise, saturating into a Lorentzian near the carrier.
* ====================================================================== */
static struct pac_harm osc_hd; /* retained oscillator conversion data */
static double *osc_Vr = NULL, *osc_Vi = NULL; /* [Ntot] retained oscillator spectrum */
static double osc_f0 = 0.0;
static int osc_node = 0, osc_valid = 0;
int
HBOSCanalyze(CKTcircuit *ckt, int oscNode, int K, int Pin, double f0seed,
double ampseed, int maxiter, double tol, int verbose)
{
int N = SMPmatSize(ckt->CKTmatrix);
int P = Pin > 0 ? Pin : ((8*K < 32) ? 32 : 8*K);
int Ntot = (2*K+1)*N, Naug = Ntot + 1;
int iter, s, i, k, c, rc = E_ITERLIM, have_hd = 0;
double f0 = f0seed, w0 = 2.0*M_PI*f0, fnorm = 0.0;
double *Vr, *Vi, *vsamp, *IRr, *IRi, *Kr, *Ki, *Jr, *Ji, *Fr, *Fi, *Ar, *Ai, *br, *bi;
struct pac_harm hd;
if (N <= 0 || K < 1 || oscNode <= 0 || oscNode > N) {
fprintf(stderr, "hbosc: bad size or oscillator node.\n"); return E_PARMVAL;
}
if (Naug > 1600) { fprintf(stderr, "hbosc: system too large (reduce K).\n"); return E_PARMVAL; }
Vr = TMALLOC(double, Ntot); Vi = TMALLOC(double, Ntot);
IRr = TMALLOC(double, Ntot); IRi = TMALLOC(double, Ntot);
Kr = TMALLOC(double, Ntot); Ki = TMALLOC(double, Ntot);
Fr = TMALLOC(double, Ntot); Fi = TMALLOC(double, Ntot);
Jr = TMALLOC(double, (size_t)Ntot*(size_t)Ntot); Ji = TMALLOC(double, (size_t)Ntot*(size_t)Ntot);
Ar = TMALLOC(double, (size_t)Naug*(size_t)Naug); Ai = TMALLOC(double, (size_t)Naug*(size_t)Naug);
br = TMALLOC(double, Naug); bi = TMALLOC(double, Naug);
vsamp = TMALLOC(double, (size_t)N*(size_t)P);
/* seed: fundamental at the osc node (amplitude A -> |V_1| = A/2, real, phase 0) */
for (i = 0; i < Ntot; i++) { Vr[i] = 0.0; Vi[i] = 0.0; }
Vr[(size_t)(1+K)*(size_t)N + (size_t)(oscNode-1)] = ampseed * 0.5;
Vr[(size_t)(K-1)*(size_t)N + (size_t)(oscNode-1)] = ampseed * 0.5;
for (iter = 0; iter < maxiter; iter++) {
w0 = 2.0*M_PI*f0;
for (s = 0; s < P; s++)
for (i = 0; i < N; i++) {
double v = 0.0;
for (k = -K; k <= K; k++) {
double ang = 2.0*M_PI*k*s/(double)P;
v += Vr[(size_t)(k+K)*(size_t)N+(size_t)i]*cos(ang)
- Vi[(size_t)(k+K)*(size_t)N+(size_t)i]*sin(ang);
}
vsamp[(size_t)s*(size_t)N+(size_t)i] = v;
}
if (hb_extract(ckt, vsamp, N, P, K, &hd, IRr, IRi)) {
fprintf(stderr, "hbosc: device extraction failed.\n"); rc = E_PARMVAL; break;
}
have_hd = 1;
pac_build_matrix(&hd, f0, 0.0, Jr, Ji);
{ /* I_C = (J - Jg) V */
struct pac_harm hg = hd;
double *Jgr = TMALLOC(double, (size_t)Ntot*(size_t)Ntot);
double *Jgi = TMALLOC(double, (size_t)Ntot*(size_t)Ntot);
hg.Cmr = TMALLOC(double, (size_t)hd.nnz*(size_t)(hd.H+1));
hg.Cmi = TMALLOC(double, (size_t)hd.nnz*(size_t)(hd.H+1));
pac_build_matrix(&hg, f0, 0.0, Jgr, Jgi);
FREE(hg.Cmr); FREE(hg.Cmi);
for (i = 0; i < Ntot; i++) {
double cr = 0, ci = 0;
for (c = 0; c < Ntot; c++) {
double ar = Jr[(size_t)i*(size_t)Ntot+(size_t)c] - Jgr[(size_t)i*(size_t)Ntot+(size_t)c];
double ai = Ji[(size_t)i*(size_t)Ntot+(size_t)c] - Jgi[(size_t)i*(size_t)Ntot+(size_t)c];
cr += ar*Vr[c] - ai*Vi[c]; ci += ar*Vi[c] + ai*Vr[c];
}
Kr[i] = cr; Ki[i] = ci;
}
FREE(Jgr); FREE(Jgi);
}
fnorm = 0.0;
for (i = 0; i < Ntot; i++) {
Fr[i] = IRr[i] + Kr[i]; Fi[i] = IRi[i] + Ki[i];
fnorm += Fr[i]*Fr[i] + Fi[i]*Fi[i];
}
fnorm = sqrt(fnorm);
if (verbose)
fprintf(stderr, "hbosc iter %2d: |F| = %.6e f0 = %.8g Hz\n", iter, fnorm, f0);
if (isnan(fnorm) || fnorm > 1e300) { pac_free_harmonics(&hd); have_hd = 0; break; }
if (fnorm < tol) { rc = OK; break; } /* converged: keep hd for retention */
/* bordered Newton: [ H d ; u*^T 0 ] [dV; dw0] = [-F; 0] */
memset(Ar, 0, (size_t)Naug*(size_t)Naug*sizeof(double));
memset(Ai, 0, (size_t)Naug*(size_t)Naug*sizeof(double));
for (i = 0; i < Ntot; i++)
for (c = 0; c < Ntot; c++) {
Ar[(size_t)i*(size_t)Naug+(size_t)c] = Jr[(size_t)i*(size_t)Ntot+(size_t)c];
Ai[(size_t)i*(size_t)Naug+(size_t)c] = Ji[(size_t)i*(size_t)Ntot+(size_t)c];
}
for (i = 0; i < Ntot; i++) { /* d = dF/dw0 = I_C / w0 */
Ar[(size_t)i*(size_t)Naug+(size_t)Ntot] = Kr[i]/w0;
Ai[(size_t)i*(size_t)Naug+(size_t)Ntot] = Ki[i]/w0;
}
for (i = 0; i < Ntot; i++) { /* row = conj(u), u_k = jk V_k */
int kk = (i/N) - K;
Ar[(size_t)Ntot*(size_t)Naug+(size_t)i] = -(double)kk * Vi[i]; /* Re conj(u) */
Ai[(size_t)Ntot*(size_t)Naug+(size_t)i] = -(double)kk * Vr[i]; /* Im conj(u) */
}
for (i = 0; i < Ntot; i++) { br[i] = -Fr[i]; bi[i] = -Fi[i]; }
br[Ntot] = 0.0; bi[Ntot] = 0.0;
pac_free_harmonics(&hd); have_hd = 0;
if (pss_csolve(Naug, Ar, Ai, br, bi)) {
fprintf(stderr, "hbosc: singular augmented system.\n"); rc = E_SINGULAR; break;
}
for (i = 0; i < Ntot; i++) { Vr[i] += br[i]; Vi[i] += bi[i]; }
f0 += br[Ntot] / (2.0*M_PI); /* dw0 = Re(x[Ntot]) */
}
if (rc == OK && have_hd) {
int numNames, error;
IFuid *nameList = NULL;
/* retain the operating point for phasenoise */
if (osc_valid) { pac_free_harmonics(&osc_hd); FREE(osc_Vr); FREE(osc_Vi); }
osc_hd = hd; /* transfer ownership of hd's arrays */
osc_Vr = TMALLOC(double, Ntot); osc_Vi = TMALLOC(double, Ntot);
memcpy(osc_Vr, Vr, (size_t)Ntot*sizeof(double));
memcpy(osc_Vi, Vi, (size_t)Ntot*sizeof(double));
osc_f0 = f0; osc_node = oscNode; osc_valid = 1;
error = CKTnames(ckt, &numNames, &nameList);
fprintf(stdout, "\nHBOSC: autonomous oscillator steady state\n"
" oscillation frequency f0 = %.9g Hz (converged, |F| = %.3e)\n"
" node harmonic frequency [Hz] |V| phase [deg]\n",
f0, fnorm);
for (i = 0; i < N; i++) {
const char *nm = (!error && i < numNames) ? (const char *) nameList[i] : "?";
for (k = 0; k <= K; k++) {
double sc = (k == 0) ? 1.0 : 2.0;
double vr = sc*Vr[(size_t)(k+K)*(size_t)N+(size_t)i];
double vi = sc*Vi[(size_t)(k+K)*(size_t)N+(size_t)i];
fprintf(stdout, " %-8s %6d %16.6e %14.6e %10.3f\n",
nm, k, k*f0, hypot(vr, vi), (k==0)?0.0:atan2(vi, vr)*180.0/M_PI);
}
}
if (nameList) tfree(nameList);
} else {
if (have_hd) pac_free_harmonics(&hd);
if (rc != OK)
fprintf(stderr, "hbosc: did not converge to an oscillation (try a better fguess/tstab).\n");
}
FREE(Vr); FREE(Vi); FREE(IRr); FREE(IRi); FREE(Kr); FREE(Ki);
FREE(Fr); FREE(Fi); FREE(Jr); FREE(Ji); FREE(Ar); FREE(Ai); FREE(br); FREE(bi); FREE(vsamp);
return rc;
}
int
PhaseNoiseAnalyze(CKTcircuit *ckt, double fstart, double fstop, int npts, int verbose)
{
struct pac_harm *hd = &osc_hd;
int N, M, Ntot, i, j, k, rc = OK, onode, ni, mi, ei;
double f0 = osc_f0;
double *Psr, *Psi, Pcar, freq, mult;
NOISEAN nj; Ndata data; JOB *oldJob;
if (!osc_valid) {
fprintf(stderr, "phasenoise: no oscillator operating point -- run `hbosc` first.\n");
return E_NOTFOUND;
}
N = hd->N; M = hd->M; Ntot = hd->Ntot; onode = osc_node;
Psr = TMALLOC(double, Ntot); Psi = TMALLOC(double, Ntot);
/* carrier power at the oscillator node: single-sided amplitude A = 2|V_1|, mean
* square A^2/2 = 2|V_1|^2 (the reference the sideband noise is measured against). */
{
size_t c1 = (size_t)(M+1)*(size_t)N + (size_t)(onode-1);
Pcar = 2.0 * (osc_Vr[c1]*osc_Vr[c1] + osc_Vi[c1]*osc_Vi[c1]);
if (Pcar <= 0.0) Pcar = 1.0;
}
/* bias the devices at the oscillator op-point (phase-0 sample) for their noise PSD */
for (j = 1; j <= N; j++) {
double v = 0.0;
for (k = -M; k <= M; k++)
v += osc_Vr[(size_t)(k+M)*(size_t)N+(size_t)(j-1)];
ckt->CKTrhsOld[j] = v;
}
ckt->CKTrhsOld[0] = 0.0;
ckt->CKTmode = (ckt->CKTmode & MODEUIC) | MODEDCOP | MODEINITSMSIG;
CKTload(ckt);
memset(&nj, 0, sizeof(nj));
nj.NstartFreq = fstart; nj.NstopFreq = fstop; nj.NnumSteps = npts;
nj.NstpType = 0; nj.NStpsSm = 0; nj.JOBname = "phasenoise";
memset(&data, 0, sizeof(data)); data.prtSummary = FALSE;
oldJob = ckt->CKTcurJob; ckt->CKTcurJob = (JOB *) &nj;
for (i = 0; i < DEVmaxnum; i++)
if (DEVices[i] && DEVices[i]->DEVnoise && ckt->CKThead[i]) {
double dummy = 0.0;
DEVices[i]->DEVnoise(N_DENS, N_OPEN, ckt->CKThead[i], ckt, &data, &dummy);
}
fprintf(stdout, "\nPhaseNoise: oscillator phase noise (f0 = %.9g Hz, carrier power %.4e)\n"
" offset [Hz] L(df) [dBc/Hz]\n", f0, Pcar);
mult = (npts > 1) ? pow(fstop/fstart, 1.0/(double)(npts-1)) : 1.0;
for (freq = fstart, i = 0; i < npts; i++, freq *= mult) {
/* adjoint of H at OFFSET f_in = df, unit at the CARRIER sideband (m=1) of the
* osc node: Psi is the transimpedance from a noise injection at every (node,
* sideband) to the carrier. As df -> 0, H(df) -> the singular limit-cycle matrix,
* so Psi (through the phase mode) blows up as 1/df -> Sv ~ 1/df^2 = phase noise. */
double *Ar = TMALLOC(double, (size_t)Ntot*(size_t)Ntot);
double *Ai = TMALLOC(double, (size_t)Ntot*(size_t)Ntot);
double Sv = 0.0;
memset(Ar, 0, (size_t)Ntot*(size_t)Ntot*sizeof(double));
memset(Ai, 0, (size_t)Ntot*(size_t)Ntot*sizeof(double));
for (ni = 0; ni <= 2*M; ni++)
for (mi = 0; mi <= 2*M; mi++) {
int dm = (ni - M) - (mi - M);
double omega = 2.0*M_PI*(freq + (double)(mi - M)*f0);
for (ei = 0; ei < hd->nnz; ei++) {
size_t hi = (size_t)ei*(size_t)(hd->H+1) + (size_t)abs(dm);
double gr = hd->Gmr[hi], gi = hd->Gmi[hi], cr = hd->Cmr[hi], ci = hd->Cmi[hi];
double er, eii; size_t row, col;
if (dm < 0) { gi = -gi; ci = -ci; }
er = gr - omega*ci; eii = gi + omega*cr;
row = (size_t)ni*(size_t)N + (size_t)(hd->rr[ei]-1);
col = (size_t)mi*(size_t)N + (size_t)(hd->cc[ei]-1);
Ar[col*(size_t)Ntot+row] += er; /* transpose: [col][row] = H^T */
Ai[col*(size_t)Ntot+row] += eii;
}
}
memset(Psr, 0, (size_t)Ntot*sizeof(double));
memset(Psi, 0, (size_t)Ntot*sizeof(double));
Psr[(size_t)(M+1)*(size_t)N + (size_t)(onode-1)] = 1.0; /* carrier sideband m=1 */
if (pss_csolve(Ntot, Ar, Ai, Psr, Psi) == 0) {
data.freq = freq; data.delFreq = 0.0; data.prtSummary = FALSE;
for (k = -M; k <= M; k++) {
double dens = 0.0;
size_t blk = (size_t)(k+M)*(size_t)N;
for (j = 1; j <= N; j++) { ckt->CKTrhs[j] = Psr[blk+(size_t)(j-1)]; ckt->CKTirhs[j] = Psi[blk+(size_t)(j-1)]; }
ckt->CKTrhs[0] = 0.0; ckt->CKTirhs[0] = 0.0;
for (ei = 0; ei < DEVmaxnum; ei++)
if (DEVices[ei] && DEVices[ei]->DEVnoise && ckt->CKThead[ei])
DEVices[ei]->DEVnoise(N_DENS, N_CALC, ckt->CKThead[ei], ckt, &data, &dens);
Sv += dens;
}
}
FREE(Ar); FREE(Ai);
fprintf(stdout, " %14.6e %12.4f\n", freq, 10.0*log10(Sv / Pcar));
if (npts == 1) break;
}
ckt->CKTcurJob = oldJob;
(void) verbose;
FREE(Psr); FREE(Psi);
return rc;
}
int
DCpss(CKTcircuit *ckt,
int restart) /* forced restart flag */

2
visualc/vngspice.vcxproj
View File

@ -887,6 +887,7 @@
<ClInclude Include="..\src\frontend\com_gnuplot.h" />
<ClInclude Include="..\src\frontend\com_hardcopy.h" />
<ClInclude Include="..\src\frontend\com_hb.h" />
<ClInclude Include="..\src\frontend\com_hbosc.h" />
<ClInclude Include="..\src\frontend\com_help.h" />
<ClInclude Include="..\src\frontend\com_history.h" />
<ClInclude Include="..\src\frontend\com_let.h" />
@ -1507,6 +1508,7 @@
<ClCompile Include="..\src\frontend\com_gnuplot.c" />
<ClCompile Include="..\src\frontend\com_hardcopy.c" />
<ClCompile Include="..\src\frontend\com_hb.c" />
<ClCompile Include="..\src\frontend\com_hbosc.c" />
<ClCompile Include="..\src\frontend\com_help.c" />
<ClCompile Include="..\src\frontend\com_history.c" />
<ClCompile Include="..\src\frontend\com_let.c" />