Fix Bug #629 - "XSPICE d_osc failures". The old code has been completely

replaced by a new design that is faster, more reliable and does not
usually insert analog breakpoints.
This commit is contained in:
Giles Atkinson 2023-07-29 07:44:22 +01:00 committed by Holger Vogt
parent 6b0ab44f1e
commit 4df2e69009
2 changed files with 185 additions and 460 deletions

View File

@ -1,479 +1,200 @@
/*.......1.........2.........3.........4.........5.........6.........7.........8
================================================================================
/* XSPICE code model for the Controlled Digital Oscillator.
* This is a complete redesign of the original version by the
* Georgia Tech team, as a fix for ngspice Bug #629 - "XSPICE d_osc failures".
*/
FILE d_osc/cfunc.mod
Public Domain
Georgia Tech Research Corporation
Atlanta, Georgia 30332
PROJECT A-8503-405
AUTHORS
24 Jul 1991 Jeffrey P. Murray
MODIFICATIONS
23 Aug 1991 Jeffrey P. Murray
30 Sep 1991 Jeffrey P. Murray
09 Nov 2022 Holger Vogt
05 Jan 2023 Robert Turnbull
SUMMARY
This file contains the model-specific routines used to
functionally describe the d_osc code model.
INTERFACES
FILE ROUTINE CALLED
CMmacros.h cm_message_send();
CM.c void *cm_analog_alloc()
void *cm_analog_get_ptr()
CMevt.c void cm_event_queue()
REFERENCED FILES
Inputs from and outputs to ARGS structure.
NON-STANDARD FEATURES
NONE
===============================================================================*/
/*=== INCLUDE FILES ====================*/
#include "d_osc.h" /* ...contains macros & type defns.
for this model. 7/24/91 - JPM */
#include <stdlib.h>
#define FACTOR 0.75 // Controls timing of next scheduled call. */
/*=== CONSTANTS ========================*/
/* PWL table entry. */
struct pwl {
double ctl, freq;
};
/* Called at end to free memory. */
/*=== MACROS ===========================*/
/*=== LOCAL VARIABLES & TYPEDEFS =======*/
typedef struct {
double *x;
double *y;
} Local_Data_t;
/*=== FUNCTION PROTOTYPE DEFINITIONS ===*/
/*==============================================================================
FUNCTION cm_d_osc()
AUTHORS
24 Jul 1991 Jeffrey P. Murray
MODIFICATIONS
30 Sep 1991 Jeffrey P. Murray
SUMMARY
This function implements the d_osc code model.
INTERFACES
FILE ROUTINE CALLED
CMmacros.h cm_message_send();
CM.c void *cm_analog_alloc()
void *cm_analog_get_ptr()
CMevt.c void cm_event_queue()
RETURNED VALUE
Returns inputs and outputs via ARGS structure.
GLOBAL VARIABLES
NONE
NON-STANDARD FEATURES
NONE
==============================================================================*/
static void cm_d_osc_callback(ARGS,
Mif_Callback_Reason_t reason)
static void callback(ARGS, Mif_Callback_Reason_t reason)
{
switch (reason) {
case MIF_CB_DESTROY: {
Local_Data_t *loc = STATIC_VAR(locdata);
if (loc) {
if (loc->x)
free(loc->x);
if(loc->y)
free(loc->y);
free(loc);
STATIC_VAR(locdata) = loc = NULL;
}
break;
} /* end of case MIF_CB_DESTROY */
} /* end of switch over reason being called */
} /* end of function cm_d_osc_callback */
struct panel_instance *instance;
if (reason == MIF_CB_DESTROY) {
struct pwl *table = STATIC_VAR(locdata);
/*=== CM_D_OSC ROUTINE ===*/
/*************************************************************
* The following is the model for the controlled digital *
* oscillator for the ATESSE Version 2.0 system. *
* *
* Created 7/24/91 J.P.Murray *
*************************************************************/
/*************************************************************
* *
* *
* <-----duty_cycle-----> *
* I *
* I t2 t3 *
* I \______________/_____ *
* I | | *
* I | | | | *
* I | | *
* I | | | | *
* I | | *
* I | | | | *
* I-----------------*-----* - - - - - - - - - -*--------- *
* t1 t4 *
* *
* *
* t2 = t1 + rise_delay *
* t4 = t3 + fall_delay *
* *
* Note that for the digital model, unlike for the *
* analog "square" model, t1 and t3 are stored and *
* adjusted values, but t2 & t4 are implied by the *
* rise and fall delays of the model, but are otherwise *
* not stored values. JPM *
* *
*************************************************************/
void cm_d_osc(ARGS)
{
double *x, /* analog input value control array */
*y, /* frequency array */
cntl_input, /* control input value */
*phase, /* instantaneous phase of the model */
*phase_old, /* previous phase of the model */
*t1, /* pointer to t1 value */
*t3, /* pointer to t3 value */
/*time1,*/ /* variable for calculating new time1 value */
/*time3,*/ /* variable for calculating new time3 value */
freq = 0.0, /* instantaneous frequency value */
dphase, /* fractional part into cycle */
duty_cycle, /* duty_cycle value */
test_double, /* testing variable */
slope; /* slope value...used to extrapolate
freq values past endpoints. */
int i, /* generic loop counter index */
cntl_size, /* control array size */
freq_size; /* frequency array size */
Local_Data_t *loc; /* Pointer to local static data, not to be included
in the state vector (save memory!) */
/**** Retrieve frequently used parameters... ****/
cntl_size = PARAM_SIZE(cntl_array);
freq_size = PARAM_SIZE(freq_array);
duty_cycle = PARAM(duty_cycle);
/* check and make sure that the control array is the
same size as the frequency array */
if(cntl_size != freq_size){
cm_message_send(d_osc_array_error);
return;
}
if (INIT) { /*** Test for INIT == TRUE. If so, allocate storage, etc. ***/
/* Allocate storage for internal variables */
cm_analog_alloc(0, sizeof(double));
cm_analog_alloc(1, sizeof(double));
cm_analog_alloc(2, sizeof(double));
/* assign internal variables */
phase = phase_old = (double *) cm_analog_get_ptr(0,0);
t1 = (double *) cm_analog_get_ptr(1,0);
t3 = (double *) cm_analog_get_ptr(2,0);
/*** allocate static storage for *loc ***/
STATIC_VAR (locdata) = calloc (1 , sizeof ( Local_Data_t ));
loc = STATIC_VAR (locdata);
CALLBACK = cm_d_osc_callback;
x = loc->x = (double *) calloc((size_t) cntl_size, sizeof(double));
if (!x) {
cm_message_send(d_osc_allocation_error);
return;
}
y = loc->y = (double *) calloc((size_t) cntl_size, sizeof(double));
if (!y) {
cm_message_send(d_osc_allocation_error);
if(x)
free(x);
return;
}
/* Retrieve x and y values. */
for (i=0; i<cntl_size; i++) {
x[i] = PARAM(cntl_array[i]);
y[i] = PARAM(freq_array[i]);
}
}
else { /*** This is not an initialization pass...retrieve storage
addresses and calculate new outputs, if required. ***/
/** Retrieve previous values... **/
/* assign internal variables */
phase = (double *) cm_analog_get_ptr(0,0);
phase_old = (double *) cm_analog_get_ptr(0,1);
t1 = (double *) cm_analog_get_ptr(1,0);
t3 = (double *) cm_analog_get_ptr(2,0);
}
switch (CALL_TYPE) {
case ANALOG: /** analog call **/
test_double = TIME;
if ( AC == ANALYSIS ) { /* this model does not function
in AC analysis mode. */
return;
}
else {
if ( 0.0 == TIME ) { /* DC analysis */
/* retrieve & normalize phase value */
*phase = PARAM(init_phase);
if ( 0 > *phase ) {
*phase = *phase + 360.0;
}
*phase = *phase / 360.0;
/* set phase value to init_phase */
*phase_old = *phase;
/* preset time values to harmless values... */
*t1 = -1;
*t3 = -1;
}
loc = STATIC_VAR (locdata);
x = loc->x;
y = loc->y;
/* Retrieve cntl_input value. */
cntl_input = INPUT(cntl_in);
/* Determine segment boundaries within which cntl_input resides */
/*** cntl_input below lowest cntl_voltage ***/
if (cntl_input <= x[0]) {
slope = (y[1] - y[0])/(x[1] - x[0]);
freq = y[0] + (cntl_input - x[0]) * slope;
}
else
/*** cntl_input above highest cntl_voltage ***/
if (cntl_input >= x[cntl_size-1]) {
slope = (y[cntl_size-1] - y[cntl_size-2]) /
(x[cntl_size-1] - x[cntl_size-2]);
freq = y[cntl_size-1] + (cntl_input - x[cntl_size-1]) * slope;
}
else { /*** cntl_input within bounds of end midpoints...
must determine position progressively & then
calculate required output. ***/
for (i=0; i<cntl_size-1; i++) {
if ( (cntl_input < x[i+1]) && (cntl_input >= x[i]) ) {
/* Interpolate to the correct frequency value */
freq = ( (cntl_input - x[i]) / (x[i+1] - x[i]) ) *
( y[i+1]-y[i] ) + y[i];
}
}
}
/*** If freq < 0.0, clamp to 1e-16 & issue a warning ***/
if ( 0.0 > freq ) {
freq = 1.0e-16;
cm_message_send(d_osc_negative_freq_error);
}
/* calculate the instantaneous phase */
*phase = *phase_old + freq * (TIME - T(1));
/* dphase is the percent into the cycle for
the period */
dphase = *phase_old - floor(*phase_old);
/* Calculate the time variables and the output value
for this iteration */
if((*t1 <= TIME) && (TIME <= *t3)) { /* output high */
*t3 = T(1) + (1 - dphase)/freq;
if(TIME < *t3) {
cm_event_queue(*t3);
}
}
else
if((*t3 <= TIME) && (TIME <= *t1)) { /* output low */
if(dphase > (1.0 - duty_cycle) ) {
dphase = dphase - 1.0;
}
*t1 = T(1) + ( (1.0 - duty_cycle) - dphase)/freq;
if(TIME < *t1) {
cm_event_queue(*t1);
}
}
else {
if(dphase > (1.0 - duty_cycle) ) {
dphase = dphase - 1.0;
}
*t1 = T(1) + ( (1.0 - duty_cycle) - dphase )/freq;
if((TIME < *t1) || (T(1) == 0)) {
cm_event_queue(*t1);
}
*t3 = T(1) + (1 - dphase)/freq;
}
}
break;
case EVENT: /** discrete call...lots to do **/
test_double = TIME;
if ( 0.0 == TIME ) { /* DC analysis...preset values,
as appropriate.... */
/* retrieve & normalize phase value */
*phase = PARAM(init_phase);
if ( 0 > *phase ) {
*phase = *phase + 360.0;
}
*phase = *phase / 360.0;
/* set phase value to init_phase */
*phase_old = *phase;
/* preset time values to harmless values... */
*t1 = -1;
*t3 = -1;
}
/* Calculate the time variables and the output value
for this iteration */
/* Output is always set to STRONG */
OUTPUT_STRENGTH(out) = STRONG;
if( *t1 == TIME ) { /* rising edge */
OUTPUT_STATE(out) = ONE;
OUTPUT_DELAY(out) = PARAM(rise_delay);
}
else {
if ( *t3 == TIME ) { /* falling edge */
OUTPUT_STATE(out) = ZERO;
OUTPUT_DELAY(out) = PARAM(fall_delay);
}
else { /* no change in output */
if ( TIME != 0.0 ) {
OUTPUT_CHANGED(out) = FALSE;
}
if ( (*t1 < TIME) && (TIME < *t3) ) {
OUTPUT_STATE(out) = ONE;
}
else {
OUTPUT_STATE(out) = ZERO;
}
}
}
break;
if (table)
free(table);
STATIC_VAR(locdata) = NULL;
}
}
/* Get the current period. */
static double get_period(double ctl, struct pwl *table, int csize)
{
double f;
int i;
for (i = 0; i < csize; ++i) {
if (table[i].ctl > ctl)
break;
}
/* Interpolation outside input range continues slope. */
if (i > 0) {
if (i == csize)
i -= 2;
else
i--;
}
f = table[i].freq +
(ctl - table[i].ctl) * ((table[i + 1].freq - table[i].freq) /
(table[i + 1].ctl - table[i].ctl));
return 1.0 / f;
}
/* The state data. */
struct state {
double last_time; // Time of last output change.
Digital_State_t last; // Last value output.
};
/* The code-model function. */
void cm_d_osc(ARGS)
{
struct pwl *table;
struct state *state;
double ctl, delta, when;
int csize, i;
csize = PARAM_SIZE(cntl_array);
if (INIT) {
/* Validate PWL table. */
for (i = 0; i < csize - 1; ++i) {
if (PARAM(cntl_array[i]) >= PARAM(cntl_array[i + 1]))
break;
}
if (i < csize - 1 || csize != PARAM_SIZE(freq_array)) {
cm_message_send("Badly-formed control table");
STATIC_VAR(locdata) = NULL;
return;
}
/* Allocate PWL table. */
table = malloc(csize * sizeof (struct pwl));
STATIC_VAR(locdata) = table;
if (!table)
return;
for (i = 0; i < csize; ++i) {
table[i].ctl = PARAM(cntl_array[i]);
table[i].freq = PARAM(freq_array[i]);
if (table[i].freq <= 0) {
cm_message_printf("Error: frequency %g is not positve, "
"value replaced by 1e-16.",
table[i].freq);
table[i].freq = 1.0e-16;
}
}
/* Allocate state data. */
cm_event_alloc(0, sizeof (struct state));
return;
}
table = STATIC_VAR(locdata);
if (!table)
return;
state = (struct state *)cm_event_get_ptr(0, 0);
if (CALL_TYPE != EVENT) {
if (TIME == 0.0) {
double phase;
/* Set initial output and state data. */
ctl = INPUT(cntl_in);
delta = get_period(ctl, table, csize);
phase = PARAM(init_phase);
phase /= 360.0;
if (phase < 0.0)
phase += 1.0;
/* When would a hypothetical previous transition have been? */
state->last_time = delta * (1.0 - PARAM(duty_cycle) - phase);
if (state->last_time < 0.0) {
state->last = ONE;
} else {
state->last = ZERO;
state->last_time = -delta * phase;
}
}
return;
}
/* Event call; either one requested previously or just before
* a time-step is accepted.
*/
if (TIME == 0.0) {
OUTPUT_STATE(out) = state->last;
OUTPUT_STRENGTH(out) = STRONG;
return;
}
/* When is the next transition due? */
ctl = INPUT(cntl_in);
delta = get_period(ctl, table, csize);
if (state->last)
delta *= PARAM(duty_cycle);
else
delta *= (1.0 - PARAM(duty_cycle));
when = state->last_time + delta;
if (TIME >= when) {
/* If the frequency rose rapidly, the transition has been missed.
* Force a shorter time-step and schedule then.
*/
cm_analog_set_temp_bkpt(state->last_time + FACTOR * delta);
OUTPUT_CHANGED(out) = FALSE;
return;
}
if (TIME >= state->last_time + FACTOR * delta) {
/* TIME is reasonably close to transition time. Request output. */
state->last_time = when;
state->last ^= ONE;
OUTPUT_STATE(out) = state->last;
OUTPUT_STRENGTH(out) = STRONG;
OUTPUT_DELAY(out) = when - TIME;
/* Request a call in the next half-cycle. */
cm_event_queue(when + FACTOR * delta);
} else {
OUTPUT_CHANGED(out) = FALSE;
if (TIME < state->last_time) {
/* Output transition pending, nothing to do. */
return;
} else {
/* Request a call nearer to transition time. */
cm_event_queue(state->last_time + FACTOR * delta);
}
}
}

View File

@ -62,6 +62,10 @@ Vector_Bounds: - -
Null_Allowed: yes yes
/* Rise and fall delay parameter are unused, but retained,
* to be compatible with an earlier version.
*/
PARAMETER_TABLE:
Parameter_Name: rise_delay fall_delay