summaryrefslogtreecommitdiffstats
path: root/src/electronics/simulation
diff options
context:
space:
mode:
authortpearson <tpearson@283d02a7-25f6-0310-bc7c-ecb5cbfe19da>2010-02-24 01:49:02 +0000
committertpearson <tpearson@283d02a7-25f6-0310-bc7c-ecb5cbfe19da>2010-02-24 01:49:02 +0000
commit5de3dd4762ca33a0f92e79ffa4fe2ff67069d531 (patch)
treebad482b7afa4cdf47422d60a5dd2c61c7e333b09 /src/electronics/simulation
downloadktechlab-5de3dd4762ca33a0f92e79ffa4fe2ff67069d531.tar.gz
ktechlab-5de3dd4762ca33a0f92e79ffa4fe2ff67069d531.zip
Added KDE3 version of ktechlab
git-svn-id: svn://anonsvn.kde.org/home/kde/branches/trinity/applications/ktechlab@1095338 283d02a7-25f6-0310-bc7c-ecb5cbfe19da
Diffstat (limited to 'src/electronics/simulation')
-rw-r--r--src/electronics/simulation/Makefile.am11
-rw-r--r--src/electronics/simulation/bjt.cpp257
-rw-r--r--src/electronics/simulation/bjt.h75
-rw-r--r--src/electronics/simulation/capacitance.cpp117
-rw-r--r--src/electronics/simulation/capacitance.h55
-rw-r--r--src/electronics/simulation/cccs.cpp89
-rw-r--r--src/electronics/simulation/cccs.h41
-rw-r--r--src/electronics/simulation/ccvs.cpp91
-rw-r--r--src/electronics/simulation/ccvs.h41
-rw-r--r--src/electronics/simulation/circuit.cpp550
-rw-r--r--src/electronics/simulation/circuit.h137
-rw-r--r--src/electronics/simulation/currentsignal.cpp67
-rw-r--r--src/electronics/simulation/currentsignal.h43
-rw-r--r--src/electronics/simulation/currentsource.cpp64
-rw-r--r--src/electronics/simulation/currentsource.h39
-rw-r--r--src/electronics/simulation/diode.cpp198
-rw-r--r--src/electronics/simulation/diode.h73
-rw-r--r--src/electronics/simulation/element.cpp193
-rw-r--r--src/electronics/simulation/element.h255
-rw-r--r--src/electronics/simulation/elementset.cpp253
-rw-r--r--src/electronics/simulation/elementset.h131
-rw-r--r--src/electronics/simulation/elementsignal.cpp64
-rw-r--r--src/electronics/simulation/elementsignal.h45
-rw-r--r--src/electronics/simulation/inductance.cpp126
-rw-r--r--src/electronics/simulation/inductance.h55
-rw-r--r--src/electronics/simulation/logic.cpp319
-rw-r--r--src/electronics/simulation/logic.h211
-rw-r--r--src/electronics/simulation/matrix.cpp546
-rw-r--r--src/electronics/simulation/matrix.h248
-rw-r--r--src/electronics/simulation/nonlinear.cpp97
-rw-r--r--src/electronics/simulation/nonlinear.h48
-rw-r--r--src/electronics/simulation/opamp.cpp77
-rw-r--r--src/electronics/simulation/opamp.h36
-rw-r--r--src/electronics/simulation/reactive.cpp36
-rw-r--r--src/electronics/simulation/reactive.h42
-rw-r--r--src/electronics/simulation/resistance.cpp95
-rw-r--r--src/electronics/simulation/resistance.h43
-rw-r--r--src/electronics/simulation/vccs.cpp93
-rw-r--r--src/electronics/simulation/vccs.h40
-rw-r--r--src/electronics/simulation/vcvs.cpp93
-rw-r--r--src/electronics/simulation/vcvs.h40
-rw-r--r--src/electronics/simulation/vec.cpp166
-rw-r--r--src/electronics/simulation/vec.h109
-rw-r--r--src/electronics/simulation/voltagepoint.cpp73
-rw-r--r--src/electronics/simulation/voltagepoint.h39
-rw-r--r--src/electronics/simulation/voltagesignal.cpp79
-rw-r--r--src/electronics/simulation/voltagesignal.h41
-rw-r--r--src/electronics/simulation/voltagesource.cpp77
-rw-r--r--src/electronics/simulation/voltagesource.h38
49 files changed, 5756 insertions, 0 deletions
diff --git a/src/electronics/simulation/Makefile.am b/src/electronics/simulation/Makefile.am
new file mode 100644
index 0000000..c45c6a0
--- /dev/null
+++ b/src/electronics/simulation/Makefile.am
@@ -0,0 +1,11 @@
+INCLUDES = -I$(top_srcdir)/src -I$(top_srcdir)/src/electronics $(all_includes)
+METASOURCES = AUTO
+noinst_LTLIBRARIES = libelements.la
+libelements_la_SOURCES = cccs.cpp ccvs.cpp circuit.cpp currentsource.cpp \
+ diode.cpp element.cpp elementset.cpp logic.cpp matrix.cpp vccs.cpp vcvs.cpp \
+ voltagesource.cpp capacitance.cpp resistance.cpp currentsignal.cpp voltagepoint.cpp \
+ voltagesignal.cpp elementsignal.cpp nonlinear.cpp reactive.cpp vec.cpp bjt.cpp opamp.cpp \
+ inductance.cpp
+noinst_HEADERS = cccs.h ccvs.h circuit.h currentsource.h diode.h element.h \
+ elementset.h logic.h matrix.h vccs.h vcvs.h voltagesource.h capacitance.h \
+ resistance.h elementsignal.h nonlinear.h reactive.h vec.h bjt.h opamp.h inductance.h
diff --git a/src/electronics/simulation/bjt.cpp b/src/electronics/simulation/bjt.cpp
new file mode 100644
index 0000000..8daab40
--- /dev/null
+++ b/src/electronics/simulation/bjt.cpp
@@ -0,0 +1,257 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "bjt.h"
+#include "diode.h"
+#include "elementset.h"
+#include "matrix.h"
+
+#include <cmath>
+using namespace std;
+
+
+//BEGIN class BJTSettings
+BJTSettings::BJTSettings()
+{
+ I_S = 1e-16;
+ N_F = 1.0;
+ N_R = 1.0;
+ B_F = 100.0;
+ B_R = 1.0;
+}
+//END class BJTSettings
+
+
+
+//BEGIN class BJTState
+BJTState::BJTState()
+{
+ reset();
+}
+
+
+void BJTState::reset()
+{
+ for ( unsigned i = 0; i < 3; ++i )
+ {
+ for ( unsigned j = 0; j < 3; ++j )
+ A[i][j] = 0.0;
+
+ I[i] = 0.0;
+ }
+}
+
+
+BJTState BJTState::operator-( const BJTState & s ) const
+{
+ BJTState newState( *this );
+
+ for ( unsigned i = 0; i < 3; ++i )
+ {
+ for ( unsigned j = 0; j < 3; ++j )
+ newState.A[i][j] -= s.A[i][j];
+
+ newState.I[i] -= s.I[i];
+ }
+
+ return newState;
+}
+//END class BJTState
+
+
+
+//BEGIN class BJT
+BJT::BJT( const bool isNPN )
+ : NonLinear()
+{
+ V_BE_prev = 0.0;
+ V_BC_prev = 0.0;
+ m_pol = isNPN ? 1 : -1;
+ m_numCNodes = 3;
+}
+
+
+BJT::~BJT()
+{
+}
+
+
+void BJT::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[0]->n(), Map::et_unstable, false );
+ }
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse( p_cnode[1]->n(), p_cnode[1]->n(), Map::et_unstable, false );
+ }
+ if ( !p_cnode[2]->isGround )
+ {
+ p_A->setUse( p_cnode[2]->n(), p_cnode[2]->n(), Map::et_unstable, false );
+ }
+
+ if ( !p_cnode[0]->isGround && !p_cnode[2]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[2]->n(), Map::et_unstable, false );
+ p_A->setUse( p_cnode[2]->n(), p_cnode[0]->n(), Map::et_unstable, false );
+ }
+ if ( !p_cnode[1]->isGround && !p_cnode[2]->isGround )
+ {
+ p_A->setUse( p_cnode[2]->n(), p_cnode[1]->n(), Map::et_unstable, false );
+ p_A->setUse( p_cnode[1]->n(), p_cnode[2]->n(), Map::et_unstable, false );
+ }
+}
+
+
+void BJT::add_initial_dc()
+{
+ V_BE_prev = 0.0;
+ V_BC_prev = 0.0;
+ m_os.reset();
+ update_dc();
+}
+
+
+void BJT::updateCurrents()
+{
+ if (!b_status)
+ return;
+
+ double V_B = p_cnode[0]->v;
+ double V_C = p_cnode[1]->v;
+ double V_E = p_cnode[2]->v;
+
+ double V_BE = (V_B - V_E) * m_pol;
+ double V_BC = (V_B - V_C) * m_pol;
+
+ double I_BE, I_BC, I_T, g_BE, g_BC, g_IF, g_IR;
+ calcIg( V_BE, V_BC, & I_BE, & I_BC, & I_T, & g_BE, & g_BC, & g_IF, & g_IR );
+
+ m_cnodeI[1] = I_BC - I_T;
+ m_cnodeI[2] = I_BE + I_T;
+ m_cnodeI[0] = -(m_cnodeI[1] + m_cnodeI[2]);
+}
+
+
+void BJT::update_dc()
+{
+ if (!b_status)
+ return;
+
+ calc_eq();
+
+ BJTState diff = m_ns - m_os;
+ for ( unsigned i = 0; i < 3; ++i )
+ {
+ for ( unsigned j = 0 ; j < 3; ++j )
+ A_g( i, j ) += diff.A[i][j];
+
+ b_i( i ) += diff.I[i];
+ }
+
+ m_os = m_ns;
+}
+
+
+void BJT::calc_eq()
+{
+ double V_B = p_cnode[0]->v;
+ double V_C = p_cnode[1]->v;
+ double V_E = p_cnode[2]->v;
+
+ double V_BE = (V_B - V_E) * m_pol;
+ double V_BC = (V_B - V_C) * m_pol;
+
+ double I_S = m_bjtSettings.I_S;
+ double N_F = m_bjtSettings.N_F;
+ double N_R = m_bjtSettings.N_R;
+
+ // adjust voltage to help convergence
+ double V_BEcrit = diodeCriticalVoltage( I_S, N_F * V_T );
+ double V_BCcrit = diodeCriticalVoltage( I_S, N_R * V_T );
+ V_BE_prev = V_BE = diodeVoltage( V_BE, V_BE_prev, V_T * N_F, V_BEcrit );
+ V_BC_prev = V_BC = diodeVoltage( V_BC, V_BC_prev, V_T * N_R, V_BCcrit );
+
+ double I_BE, I_BC, I_T, g_BE, g_BC, g_IF, g_IR;
+ calcIg( V_BE, V_BC, & I_BE, & I_BC, & I_T, & g_BE, & g_BC, & g_IF, & g_IR );
+
+ double I_eq_B = I_BE - V_BE * g_BE;
+ double I_eq_C = I_BC - V_BC * g_BC;
+ double I_eq_E = I_T - V_BE * g_IF + V_BC * g_IR;
+
+ m_ns.A[0][0] = g_BC + g_BE;
+ m_ns.A[0][1] = -g_BC;
+ m_ns.A[0][2] = -g_BE;
+
+ m_ns.A[1][0] = -g_BC + (g_IF - g_IR);
+ m_ns.A[1][1] = g_IR + g_BC;
+ m_ns.A[1][2] = -g_IF;
+
+ m_ns.A[2][0] = -g_BE - (g_IF - g_IR);
+ m_ns.A[2][1] = -g_IR;
+ m_ns.A[2][2] = g_BE + g_IF;
+
+ m_ns.I[0] = (-I_eq_B - I_eq_C) * m_pol;
+ m_ns.I[1] = (+I_eq_C - I_eq_E) * m_pol;
+ m_ns.I[2] = (+I_eq_B + I_eq_E) * m_pol;
+}
+
+
+void BJT::calcIg( double V_BE, double V_BC, double * I_BE, double * I_BC, double * I_T, double * g_BE, double * g_BC, double * g_IF, double * g_IR )
+{
+ double I_S = m_bjtSettings.I_S;
+ double N_F = m_bjtSettings.N_F;
+ double N_R = m_bjtSettings.N_R;
+ double B_F = m_bjtSettings.B_F;
+ double B_R = m_bjtSettings.B_R;
+
+ // base-emitter diodes
+ double g_tiny = (V_BE < (-10 * V_T * N_F)) ? I_S : 0;
+
+ double I_F;
+ diodeJunction( V_BE, I_S, V_T * N_F, & I_F, g_IF );
+
+ double I_BEI = I_F / B_F;
+ double g_BEI = *g_IF / B_F;
+ double I_BEN = g_tiny * V_BE;
+ double g_BEN = g_tiny;
+ *I_BE = I_BEI + I_BEN;
+ *g_BE = g_BEI + g_BEN;
+
+ // base-collector diodes
+ g_tiny = (V_BC < (-10 * V_T * N_R)) ? I_S : 0;
+
+ double I_R;
+ diodeJunction( V_BC, I_S, V_T * N_R, & I_R, g_IR );
+
+ double I_BCI = I_R / B_R;
+ double g_BCI = *g_IR / B_R;
+ double I_BCN = g_tiny * V_BC;
+ double g_BCN = g_tiny;
+ *I_BC = I_BCI + I_BCN;
+ *g_BC = g_BCI + g_BCN;
+
+ *I_T = I_F - I_R;
+}
+
+
+void BJT::setBJTSettings( const BJTSettings & settings )
+{
+ m_bjtSettings = settings;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+}
+//END class BJT
+
+
diff --git a/src/electronics/simulation/bjt.h b/src/electronics/simulation/bjt.h
new file mode 100644
index 0000000..52022aa
--- /dev/null
+++ b/src/electronics/simulation/bjt.h
@@ -0,0 +1,75 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef BJT_H
+#define BJT_H
+
+#include "matrix.h"
+#include "nonlinear.h"
+
+class BJTState
+{
+ public:
+ BJTState();
+ void reset();
+
+ BJTState operator-( const BJTState & s ) const;
+
+ double A[3][3];
+ double I[3];
+};
+
+
+class BJTSettings
+{
+ public:
+ BJTSettings();
+
+ double I_S; ///< saturation current
+ double N_F; ///< forward emission coefficient
+ double N_R; ///< reverse emission coefficient
+ double B_F; ///< forward beta
+ double B_R; ///< reverse beta
+};
+
+
+/**
+@author David Saxton
+*/
+class BJT : public NonLinear
+{
+ public:
+ BJT( bool isNPN );
+ virtual ~BJT();
+
+ virtual Type type() const { return Element_BJT; }
+ virtual void update_dc();
+ virtual void add_initial_dc();
+ virtual void add_map();
+ BJTSettings settings() const { return m_bjtSettings; }
+ void setBJTSettings( const BJTSettings & settings );
+
+ protected:
+ virtual void updateCurrents();
+ /**
+ * Calculates the new BJTState from the voltages on the nodes.
+ */
+ void calc_eq();
+
+ void calcIg( double V_BE, double V_BC, double * I_BE, double * I_BC, double * I_T, double * g_BE, double * g_BC, double * g_IF, double * g_IR );
+
+ BJTState m_os;
+ BJTState m_ns;
+ int m_pol;
+ double V_BE_prev, V_BC_prev;
+ BJTSettings m_bjtSettings;
+};
+
+#endif
diff --git a/src/electronics/simulation/capacitance.cpp b/src/electronics/simulation/capacitance.cpp
new file mode 100644
index 0000000..9087c7f
--- /dev/null
+++ b/src/electronics/simulation/capacitance.cpp
@@ -0,0 +1,117 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "capacitance.h"
+#include "matrix.h"
+
+Capacitance::Capacitance( const double capacitance, const double delta )
+ : Reactive(delta)
+{
+ m_cap = capacitance;
+ g_eq_old = i_eq_old = 0.;
+ m_numCNodes = 2;
+ setMethod( Capacitance::m_euler );
+}
+
+Capacitance::~Capacitance()
+{
+}
+
+void Capacitance::setCapacitance( const double c )
+{
+ m_cap = c;
+}
+
+void Capacitance::add_initial_dc()
+{
+ // We don't need to do anything here, as time_step() will do that for us,
+ // apart from to make sure our old values are 0
+ g_eq_old = i_eq_old = 0.;
+}
+
+void Capacitance::updateCurrents()
+{
+ if (!b_status) return;
+ const double r_i = (p_cnode[0]->v-p_cnode[1]->v)*g_eq_old;
+ m_cnodeI[0] = -i_eq_old-r_i;
+ m_cnodeI[1] = -m_cnodeI[0];
+}
+
+
+void Capacitance::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[0]->n(), Map::et_unstable, false );
+ }
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse( p_cnode[1]->n(), p_cnode[1]->n(), Map::et_unstable, false );
+ }
+
+ if ( !p_cnode[0]->isGround && !p_cnode[1]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[1]->n(), Map::et_unstable, false );
+ p_A->setUse( p_cnode[1]->n(), p_cnode[0]->n(), Map::et_unstable, false );
+ }
+}
+
+
+void Capacitance::time_step()
+{
+ if (!b_status) return;
+
+ double v = p_cnode[0]->v-p_cnode[1]->v;
+ double i_eq_new = 0.0, g_eq_new = 0.0;
+
+ if ( m_method == Capacitance::m_euler )
+ {
+ g_eq_new = m_cap/m_delta;
+ i_eq_new = -v*g_eq_new;
+ }
+ else if ( m_method == Capacitance::m_trap ) {
+ // TODO Implement + test trapezoidal method
+ g_eq_new = 2.*m_cap/m_delta;
+ }
+
+ if ( g_eq_old != g_eq_new )
+ {
+ A_g( 0, 0 ) += g_eq_new-g_eq_old;
+ A_g( 1, 1 ) += g_eq_new-g_eq_old;
+ A_g( 0, 1 ) -= g_eq_new-g_eq_old;
+ A_g( 1, 0 ) -= g_eq_new-g_eq_old;
+ }
+
+ if ( i_eq_new != i_eq_old )
+ {
+ b_i( 0 ) -= i_eq_new-i_eq_old;
+ b_i( 1 ) += i_eq_new-i_eq_old;
+ }
+
+ g_eq_old = g_eq_new;
+ i_eq_old = i_eq_new;
+}
+
+bool Capacitance::updateStatus()
+{
+ b_status = Reactive::updateStatus();
+ if ( m_method == Capacitance::m_none ) b_status = false;
+ return b_status;
+}
+
+void Capacitance::setMethod( Method m )
+{
+ m_method = m;
+ updateStatus();
+}
+
+
diff --git a/src/electronics/simulation/capacitance.h b/src/electronics/simulation/capacitance.h
new file mode 100644
index 0000000..ccc083d
--- /dev/null
+++ b/src/electronics/simulation/capacitance.h
@@ -0,0 +1,55 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef CAPACITANCE_H
+#define CAPACITANCE_H
+
+#include "reactive.h"
+
+/**
+@author David Saxton
+@short Capacitance
+*/
+class Capacitance : public Reactive
+{
+public:
+ enum Method
+ {
+ m_none, // None
+ m_euler, // Backward Euler
+ m_trap // Trapezoidal (currently unimplemented)
+ };
+ Capacitance( const double capacitance, const double delta );
+ virtual ~Capacitance();
+
+ virtual Type type() const { return Element_Capacitance; }
+ /**
+ * Set the stepping use for numerical integration of capacitance,
+ * and the interval between successive updates
+ */
+ void setMethod( Method m );
+ virtual void time_step();
+ virtual void add_initial_dc();
+ void setCapacitance( const double c );
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual bool updateStatus();
+
+private:
+ double m_cap; // Capacitance
+ Method m_method; // Method of integration
+
+ double g_eq_old;
+ double i_eq_old;
+};
+
+#endif
diff --git a/src/electronics/simulation/cccs.cpp b/src/electronics/simulation/cccs.cpp
new file mode 100644
index 0000000..9725735
--- /dev/null
+++ b/src/electronics/simulation/cccs.cpp
@@ -0,0 +1,89 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementset.h"
+#include "matrix.h"
+#include "cccs.h"
+
+CCCS::CCCS( const double gain )
+ : Element::Element()
+{
+ m_g = gain;
+ m_numCBranches = 2;
+ m_numCNodes = 4;
+}
+
+CCCS::~CCCS()
+{
+}
+
+void CCCS::setGain( const double g )
+{
+ if ( m_g == g )
+ return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ m_g = g;
+ add_initial_dc();
+}
+
+
+void CCCS::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_b( p_cnode[0]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse_b( p_cnode[1]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[1]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[2]->isGround )
+ {
+ p_A->setUse_b( p_cnode[2]->n(), p_cbranch[1]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[3]->isGround )
+ {
+ p_A->setUse_b( p_cnode[3]->n(), p_cbranch[1]->n(), Map::et_constant, true );
+ }
+ p_A->setUse_d( p_cbranch[1]->n(), p_cbranch[0]->n(), Map::et_stable, true );
+ p_A->setUse_d( p_cbranch[1]->n(), p_cbranch[1]->n(), Map::et_constant, true );
+}
+
+void CCCS::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ A_b( 0, 0 ) = 1;
+ A_c( 0, 0 ) = 1;
+ A_b( 1, 0 ) = -1;
+ A_c( 0, 1 ) = -1;
+ A_b( 2, 1 ) = 1;
+ A_b( 3, 1 ) = -1;
+ A_d( 1, 0 ) = -m_g;
+ A_d( 1, 1 ) = 1;
+}
+
+void CCCS::updateCurrents()
+{
+ if (!b_status) return;
+ m_cnodeI[1] = p_cbranch[0]->i;
+ m_cnodeI[0] = -m_cnodeI[1];
+ m_cnodeI[3] = p_cbranch[1]->i;
+ m_cnodeI[2] = -m_cnodeI[3];
+}
+
diff --git a/src/electronics/simulation/cccs.h b/src/electronics/simulation/cccs.h
new file mode 100644
index 0000000..ba86e9e
--- /dev/null
+++ b/src/electronics/simulation/cccs.h
@@ -0,0 +1,41 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef CCCS_H
+#define CCCS_H
+
+#include "element.h"
+
+/**
+CNodes n0 and n1 are used for the current control.
+CNodes n2 and n3 are used for the current output.
+Branches b0 and b1 are for control and output
+@short Current Controlled Current Source
+@author David Saxton
+*/
+class CCCS : public Element
+{
+public:
+ CCCS( const double gain );
+ virtual ~CCCS();
+
+ virtual Type type() const { return Element_CCCS; }
+ void setGain( const double g );
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_g; // Conductance
+};
+
+#endif
diff --git a/src/electronics/simulation/ccvs.cpp b/src/electronics/simulation/ccvs.cpp
new file mode 100644
index 0000000..fc12bf5
--- /dev/null
+++ b/src/electronics/simulation/ccvs.cpp
@@ -0,0 +1,91 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementset.h"
+#include "matrix.h"
+#include "ccvs.h"
+
+CCVS::CCVS( const double gain )
+ : Element::Element()
+{
+ m_g = gain;
+ m_numCBranches = 2;
+ m_numCNodes = 4;
+}
+
+CCVS::~CCVS()
+{
+}
+
+void CCVS::setGain( const double g )
+{
+ if ( m_g == g )
+ return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ m_g = g;
+ add_initial_dc();
+}
+
+
+void CCVS::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_b( p_cnode[0]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse_b( p_cnode[1]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[1]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[2]->isGround )
+ {
+ p_A->setUse_b( p_cnode[2]->n(), p_cbranch[1]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[1]->n(), p_cnode[2]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[3]->isGround )
+ {
+ p_A->setUse_b( p_cnode[3]->n(), p_cbranch[1]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[1]->n(), p_cnode[3]->n(), Map::et_constant, true );
+ }
+ p_A->setUse_d( p_cbranch[1]->n(), p_cbranch[0]->n(), Map::et_stable, true );
+}
+
+
+void CCVS::add_initial_dc()
+{
+ if (!b_status) return;
+
+ A_b( 0, 0 ) = 1;
+ A_c( 0, 0 ) = 1;
+ A_b( 1, 0 ) = -1;
+ A_c( 0, 1 ) = -1;
+ A_b( 2, 1 ) = 1;
+ A_c( 1, 2 ) = 1;
+ A_b( 3, 1 ) = -1;
+ A_c( 1, 3 ) = -1;
+ A_d( 1, 0 ) = -m_g;
+}
+
+void CCVS::updateCurrents()
+{
+ if (!b_status) return;
+ m_cnodeI[1] = p_cbranch[0]->i;
+ m_cnodeI[0] = -m_cnodeI[1];
+ m_cnodeI[3] = p_cbranch[0]->i;
+ m_cnodeI[2] = -m_cnodeI[3];
+}
+
diff --git a/src/electronics/simulation/ccvs.h b/src/electronics/simulation/ccvs.h
new file mode 100644
index 0000000..bcb1ac0
--- /dev/null
+++ b/src/electronics/simulation/ccvs.h
@@ -0,0 +1,41 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef CCVS_H
+#define CCVS_H
+
+#include "element.h"
+
+/**
+CNodes n0 and n1 are used for the current control.
+CNodes n2 and n3 are used for the voltage output.
+Branches b0 and b1 are for control and output
+@short Current Controlled Voltage Source
+@author David Saxton
+*/
+class CCVS : public Element
+{
+public:
+ CCVS( const double gain );
+ virtual ~CCVS();
+
+ virtual Type type() const { return Element_CCVS; }
+ void setGain( const double g );
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_g; // Conductance
+};
+
+#endif
diff --git a/src/electronics/simulation/circuit.cpp b/src/electronics/simulation/circuit.cpp
new file mode 100644
index 0000000..c152756
--- /dev/null
+++ b/src/electronics/simulation/circuit.cpp
@@ -0,0 +1,550 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include <vector>
+#include "circuit.h"
+#include "circuitdocument.h"
+#include "element.h"
+#include "elementset.h"
+#include "logic.h"
+#include "matrix.h"
+#include "nonlinear.h"
+#include "pin.h"
+#include "reactive.h"
+#include "wire.h"
+
+
+#include <cmath>
+#include <map>
+
+typedef std::multimap<int, PinList> PinListMap;
+
+//BEGIN class Circuit
+Circuit::Circuit()
+{
+ m_bCanAddChanged = true;
+ m_pNextChanged[0] = m_pNextChanged[1] = 0l;
+ m_logicOutCount = 0;
+ m_bCanCache = false;
+ m_pLogicOut = 0l;
+ m_elementSet = new ElementSet( this, 0, 0 );
+ m_cnodeCount = m_branchCount = -1;
+ m_prepNLCount = 0;
+ m_pLogicCacheBase = new LogicCacheNode;
+}
+
+Circuit::~Circuit()
+{
+ delete m_elementSet;
+ delete m_pLogicCacheBase;
+ delete[] m_pLogicOut;
+}
+
+
+void Circuit::addPin( Pin *node )
+{
+ if ( m_pinList.contains(node) ) return;
+ m_pinList.append(node);
+}
+
+void Circuit::addElement( Element *element )
+{
+ if ( m_elementList.contains(element) ) return;
+ m_elementList.append(element);
+}
+
+bool Circuit::contains( Pin *node )
+{
+ return m_pinList.contains(node);
+}
+
+
+// static function
+int Circuit::identifyGround( PinList nodeList, int *highest )
+{
+ // What this function does:
+ // We are given a list of pins. First, we divide them into groups of pins
+ // that are directly connected to each other (e.g. through wires or
+ // switches). Then, each group of connected pins is looked at to find the
+ // pin with the highest "ground priority", and this is taken to be
+ // the priority of the group. The highest ground priority from all the
+ // groups is recorded. If the highest ground priority found is the maximum,
+ // then all the pins in groups with this priority are marked as ground
+ // (their eq-id is set to -1). Otherwise, the first group of pins with the
+ // highest ground priority found is marked as ground, and all others are
+ // marked as non ground (their eq-id is set to 0).
+
+ int temp_highest;
+ if (!highest)
+ highest = &temp_highest;
+
+ // Now to give all the Pins ids
+ PinListMap eqs;
+ while ( !nodeList.isEmpty() )
+ {
+ PinList associated;
+ PinList nodes;
+ Pin *node = *nodeList.begin();
+ recursivePinAdd( node, &nodeList, &associated, &nodes );
+ if ( nodes.size() > 0 )
+ {
+ eqs.insert( std::make_pair( associated.size(), nodes ) );
+ }
+ }
+
+
+ // Now, we want to look through the associated Pins,
+ // to find the ones with the highest "Ground Priority". Anything with a lower
+ // priority than Pin::gt_never will not be considered
+ *highest = Pin::gt_never; // The highest priority found so far
+ int numGround = 0; // The number of node groups found with that priority
+ const PinListMap::iterator eqsEnd = eqs.end();
+ for ( PinListMap::iterator it = eqs.begin(); it != eqsEnd; ++it )
+ {
+ int highPri = Pin::gt_never; // The highest priority found in these group of nodes
+ const PinList::iterator send = it->second.end();
+ for ( PinList::iterator sit = it->second.begin(); sit != send; ++sit )
+ {
+ if ( (*sit)->groundType() < highPri )
+ highPri = (*sit)->groundType();
+ }
+
+ if ( highPri == *highest )
+ numGround++;
+
+ else if ( highPri < *highest )
+ {
+ numGround = 1;
+ *highest = highPri;
+ }
+ }
+
+ if ( *highest == Pin::gt_never )
+ {
+ (*highest)--;
+ numGround=0;
+ }
+ // If there are no Always Ground nodes, then we only want to set one of the nodes as ground
+ else if ( *highest > Pin::gt_always )
+ numGround = 1;
+
+
+ // Now, we can give the nodes their cnode ids, or tell them they are ground
+ bool foundGround = false; // This is only used when we don't have a Always ground node
+ for ( PinListMap::iterator it = eqs.begin(); it != eqsEnd; ++it )
+ {
+ bool ground = false;
+ const PinList::iterator send = it->second.end();
+ for ( PinList::iterator sit = it->second.begin(); sit != send; ++sit )
+ {
+ ground |= (*sit)->groundType() <= (*highest);
+ }
+ if ( ground && (!foundGround || *highest == Pin::gt_always ) )
+ {
+ for ( PinList::iterator sit = it->second.begin(); sit != send; ++sit )
+ {
+ (*sit)->setEqId(-1);
+ }
+ foundGround = true;
+ }
+ else
+ {
+ for ( PinList::iterator sit = it->second.begin(); sit != send; ++sit )
+ {
+ (*sit)->setEqId(0);
+ }
+ }
+ }
+
+ return numGround;
+}
+
+
+void Circuit::init()
+{
+ m_branchCount = 0;
+
+ const ElementList::iterator listEnd = m_elementList.end();
+ for ( ElementList::iterator it = m_elementList.begin(); it != listEnd; ++it )
+ {
+ m_branchCount += (*it)->numCBranches();
+ }
+
+ // Now to give all the Pins ids
+ int groundCount = 0;
+ PinListMap eqs;
+ PinList unassignedNodes = m_pinList;
+ while ( !unassignedNodes.isEmpty() )
+ {
+ PinList associated;
+ PinList nodes;
+ Pin *node = *unassignedNodes.begin();
+ if ( recursivePinAdd( node, &unassignedNodes, &associated, &nodes ) ) {
+ groundCount++;
+ }
+ if ( nodes.size() > 0 ) {
+ eqs.insert( std::make_pair( associated.size(), nodes ) );
+ }
+ }
+
+ m_cnodeCount = eqs.size() - groundCount;
+
+ delete m_pLogicCacheBase;
+ m_pLogicCacheBase = 0l;
+
+ delete m_elementSet;
+ m_elementSet = new ElementSet( this, m_cnodeCount, m_branchCount );
+
+ // Now, we can give the nodes their cnode ids, or tell them they are ground
+ Vector *x = m_elementSet->x();
+ int i=0;
+ const PinListMap::iterator eqsEnd = eqs.end();
+ for ( PinListMap::iterator it = eqs.begin(); it != eqsEnd; ++it )
+ {
+ bool foundGround = false;
+
+ const PinList::iterator sEnd = it->second.end();
+ for ( PinList::iterator sit = it->second.begin(); sit != sEnd; ++sit )
+ foundGround |= (*sit)->eqId() == -1;
+
+ if ( foundGround )
+ continue;
+
+ bool foundEnergyStoragePin = false;
+
+ for ( PinList::iterator sit = it->second.begin(); sit != sEnd; ++sit )
+ {
+ (*sit)->setEqId(i);
+
+ bool energyStorage = false;
+ const ElementList elements = (*sit)->elements();
+ ElementList::const_iterator elementsEnd = elements.end();
+ for ( ElementList::const_iterator it = elements.begin(); it != elementsEnd; ++it )
+ {
+ if ( !*it )
+ continue;
+
+ if ( ((*it)->type() == Element::Element_Capacitance)
+ || ((*it)->type() == Element::Element_Inductance) )
+ {
+ energyStorage = true;
+ break;
+ }
+ }
+
+ // A pin attached to an energy storage pin overrides one that doesn't.
+ // If the two pins have equal status with in this regard, we pick the
+ // one with the highest absolute voltage on it.
+
+ if ( foundEnergyStoragePin && !energyStorage )
+ continue;
+
+ double v = (*sit)->voltage();
+
+ if ( energyStorage && !foundEnergyStoragePin )
+ {
+ foundEnergyStoragePin = true;
+ (*x)[i] = v;
+ continue;
+ }
+
+ if ( std::abs(v) > std::abs( (*x)[i] ) )
+ (*x)[i] = v;
+ }
+ i++;
+ }
+
+
+ // And add the elements to the elementSet
+ for ( ElementList::iterator it = m_elementList.begin(); it != listEnd; ++it )
+ {
+ // We don't want the element to prematurely try to do anything,
+ // as it doesn't know its actual cnode ids yet
+ (*it)->setCNodes();
+ (*it)->setCBranches();
+ m_elementSet->addElement(*it);
+ }
+ // And give the branch ids to the elements
+ i=0;
+ for ( ElementList::iterator it = m_elementList.begin(); it != listEnd; ++it )
+ {
+ switch ( (*it)->numCBranches() )
+ {
+ case 0:
+ break;
+ case 1:
+ (*it)->setCBranches( i );
+ i += 1;
+ break;
+ case 2:
+ (*it)->setCBranches( i, i+1 );
+ i += 2;
+ break;
+ case 3:
+ (*it)->setCBranches( i, i+1, i+2 );
+ i += 3;
+ break;
+ default:
+ // What the?!
+ break;
+ }
+ }
+}
+
+
+void Circuit::initCache()
+{
+ m_elementSet->updateInfo();
+
+ m_bCanCache = true;
+ m_logicOutCount = 0;
+
+ delete[] m_pLogicOut;
+ m_pLogicOut = 0l;
+
+ delete m_pLogicCacheBase;
+ m_pLogicCacheBase = 0l;
+
+ const ElementList::iterator end = m_elementList.end();
+ for ( ElementList::iterator it = m_elementList.begin(); it != end && m_bCanCache; ++it )
+ {
+ switch ( (*it)->type() )
+ {
+ case Element::Element_BJT:
+ case Element::Element_CCCS:
+ case Element::Element_CCVS:
+ case Element::Element_CurrentSource:
+ case Element::Element_Diode:
+ case Element::Element_LogicIn:
+ case Element::Element_OpAmp:
+ case Element::Element_Resistance:
+ case Element::Element_VCCS:
+ case Element::Element_VCVS:
+ case Element::Element_VoltagePoint:
+ case Element::Element_VoltageSource:
+ {
+ break;
+ }
+
+ case Element::Element_LogicOut:
+ {
+ m_logicOutCount++;
+ break;
+ }
+
+ case Element::Element_CurrentSignal:
+ case Element::Element_VoltageSignal:
+ case Element::Element_Capacitance:
+ case Element::Element_Inductance:
+ {
+ m_bCanCache = false;
+ break;
+ }
+ }
+ }
+
+ if ( !m_bCanCache )
+ return;
+
+ m_pLogicOut = new LogicOut*[m_logicOutCount];
+ unsigned i = 0;
+ for ( ElementList::iterator it = m_elementList.begin(); it != end && m_bCanCache; ++it )
+ {
+ if ( (*it)->type() == Element::Element_LogicOut )
+ m_pLogicOut[i++] = static_cast<LogicOut*>(*it);
+ }
+
+ m_pLogicCacheBase = new LogicCacheNode;
+}
+
+
+void Circuit::setCacheInvalidated()
+{
+ if (m_pLogicCacheBase)
+ {
+ delete m_pLogicCacheBase->high;
+ m_pLogicCacheBase->high = 0l;
+
+ delete m_pLogicCacheBase->low;
+ m_pLogicCacheBase->low = 0l;
+
+ delete m_pLogicCacheBase->data;
+ m_pLogicCacheBase->data = 0l;
+ }
+}
+
+
+void Circuit::cacheAndUpdate()
+{
+ LogicCacheNode * node = m_pLogicCacheBase;
+ for ( unsigned i = 0; i < m_logicOutCount; i++ )
+ {
+ if ( m_pLogicOut[i]->outputState() )
+ {
+ if (!node->high)
+ node->high = new LogicCacheNode;
+
+ node = node->high;
+ }
+ else
+ {
+ if (!node->low)
+ node->low = new LogicCacheNode;
+
+ node = node->low;
+ }
+ }
+
+ if ( node->data )
+ {
+ (*m_elementSet->x()) = *node->data;
+ m_elementSet->updateInfo();
+ return;
+ }
+
+ if ( m_elementSet->containsNonLinear() )
+ m_elementSet->doNonLinear( 150, 1e-10, 1e-13 );
+ else
+ m_elementSet->doLinear(true);
+
+ node->data = new Vector( m_elementSet->x()->size() );
+ *node->data = *m_elementSet->x();
+}
+
+
+void Circuit::createMatrixMap()
+{
+ m_elementSet->createMatrixMap();
+}
+
+
+bool Circuit::recursivePinAdd( Pin *node, PinList *unassignedNodes, PinList *associated, PinList *nodes )
+{
+ if ( !unassignedNodes->contains(node) )
+ return false;
+ unassignedNodes->remove(node);
+
+ bool foundGround = node->eqId() == -1;
+
+ const PinList circuitDependentPins = node->circuitDependentPins();
+ const PinList::const_iterator dEnd = circuitDependentPins.end();
+ for ( PinList::const_iterator it = circuitDependentPins.begin(); it != dEnd; ++it )
+ {
+ if ( !associated->contains(*it) )
+ associated->append(*it);
+ }
+
+ nodes->append(node);
+
+ const PinList localConnectedPins = node->localConnectedPins();
+ const PinList::const_iterator end = localConnectedPins.end();
+ for ( PinList::const_iterator it = localConnectedPins.begin(); it != end; ++it )
+ foundGround |= recursivePinAdd( *it, unassignedNodes, associated, nodes );
+
+ return foundGround;
+}
+
+
+void Circuit::doNonLogic()
+{
+ if ( !m_elementSet || m_cnodeCount+m_branchCount <= 0 )
+ return;
+
+ if (m_bCanCache)
+ {
+ if ( !m_elementSet->b()->isChanged() && !m_elementSet->matrix()->isChanged() )
+ return;
+ cacheAndUpdate();
+ updateNodalVoltages();
+ m_elementSet->b()->setUnchanged();
+ return;
+ }
+
+ stepReactive();
+ if ( m_elementSet->containsNonLinear() )
+ {
+ m_elementSet->doNonLinear( 10, 1e-9, 1e-12 );
+ updateNodalVoltages();
+ }
+ else
+ {
+ if ( m_elementSet->doLinear(true) )
+ updateNodalVoltages();
+ }
+}
+
+
+void Circuit::stepReactive()
+{
+ ElementList::iterator listEnd = m_elementList.end();
+ for ( ElementList::iterator it = m_elementList.begin(); it != listEnd; ++it )
+ {
+ Element * const e = *it;
+ if ( e && e->isReactive() )
+ (static_cast<Reactive*>(e))->time_step();
+ }
+}
+
+
+void Circuit::updateNodalVoltages()
+{
+ CNode **_cnodes = m_elementSet->cnodes();
+
+ const PinList::iterator endIt = m_pinList.end();
+ for ( PinList::iterator it = m_pinList.begin(); it != endIt; ++it )
+ {
+ Pin * const node = *it;
+ int i = node->eqId();
+ if ( i == -1 )
+ node->setVoltage(0.);
+ else
+ {
+ const double v = _cnodes[i]->v;
+ node->setVoltage( std::isfinite(v)?v:0. );
+ }
+ }
+}
+
+
+void Circuit::updateCurrents()
+{
+ ElementList::iterator listEnd = m_elementList.end();
+ for ( ElementList::iterator it = m_elementList.begin(); it != listEnd; ++it )
+ {
+ (*it)->updateCurrents();
+ }
+}
+
+void Circuit::displayEquations()
+{
+ m_elementSet->displayEquations();
+}
+//END class Circuit
+
+
+
+//BEGIN class LogicCacheNode
+LogicCacheNode::LogicCacheNode()
+{
+ low = 0l;
+ high = 0l;
+ data = 0l;
+}
+
+
+LogicCacheNode::~LogicCacheNode()
+{
+ delete low;
+ delete high;
+ delete data;
+}
+//END class LogicCacheNode
+
+
diff --git a/src/electronics/simulation/circuit.h b/src/electronics/simulation/circuit.h
new file mode 100644
index 0000000..2455edc
--- /dev/null
+++ b/src/electronics/simulation/circuit.h
@@ -0,0 +1,137 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef CIRCUIT_H
+#define CIRCUIT_H
+
+#include <qguardedptr.h>
+#include "qstringlist.h"
+#include "qvaluelist.h"
+
+#include "elementset.h"
+
+class CircuitDocument;
+class Wire;
+class Pin;
+class Element;
+class LogicOut;
+
+typedef QValueList<QGuardedPtr<Pin> > PinList;
+typedef QValueList<Element*> ElementList;
+
+
+class LogicCacheNode
+{
+ public:
+ LogicCacheNode();
+ ~LogicCacheNode();
+
+ LogicCacheNode * high;
+ LogicCacheNode * low;
+ Vector * data;
+};
+
+
+/**
+Usage of this class (usually invoked from CircuitDocument):
+(1) Add Wires, Pins and Elements to the class as appropriate
+(2) Call init to initialize the simulation
+(3) Control the simulation with step()
+
+This class can be considered a bridge between the gui-tainted CircuitDocument - specific
+to this implementation, and the pure untainted ElementSet. Please keep it that way.
+
+@short Simulates a collection of components
+@author David Saxton
+*/
+class Circuit
+{
+ public:
+ Circuit();
+ ~Circuit();
+
+ void addPin( Pin *node );
+ void addElement( Element *element );
+
+ bool contains( Pin *node );
+ bool containsNonLinear() const { return m_elementSet->containsNonLinear(); }
+
+ void init();
+ /**
+ * Called after everything else has been setup - before doNonLogic or
+ * doLogic are called for the first time. Preps the circuit.
+ */
+ void initCache();
+ /**
+ * Marks all cached results as invalidated and removes them.
+ */
+ void setCacheInvalidated();
+ /**
+ * Solves for non-logic elements
+ */
+ void doNonLogic();
+ /**
+ * Solves for logic elements (i.e just does fbSub)
+ */
+ void doLogic() { m_elementSet->doLinear(false); }
+
+ void displayEquations();
+ void updateCurrents();
+
+ void createMatrixMap();
+ /**
+ * This will identify the ground node and non-ground nodes in the given set.
+ * Ground will be given the eqId -1, non-ground of 0.
+ * @param highest The highest ground type of the groundnodes found. If no
+ ground nodes were found, this will be (gt_never-1).
+ * @returns the number of ground nodes. If all nodes are at or below the
+ * gt_never threshold, then this will be zero.
+ */
+ static int identifyGround( PinList nodeList, int *highest = 0l );
+
+ void setNextChanged( Circuit * circuit, unsigned char chain ) { m_pNextChanged[chain] = circuit; }
+ Circuit * nextChanged( unsigned char chain ) const { return m_pNextChanged[chain]; }
+ void setCanAddChanged( bool canAdd ) { m_bCanAddChanged = canAdd; }
+ bool canAddChanged() const { return m_bCanAddChanged; }
+
+ protected:
+ void cacheAndUpdate();
+ /**
+ * Update the nodal voltages from those calculated in ElementSet
+ */
+ void updateNodalVoltages();
+ /**
+ * Step the reactive elements.
+ */
+ void stepReactive();
+ /**
+ * Returns true if any of the nodes are ground
+ */
+ static bool recursivePinAdd( Pin *node, PinList *unassignedNodes, PinList *associated, PinList *nodes );
+
+ int m_cnodeCount;
+ int m_branchCount;
+ int m_prepNLCount; // Count until next m_elementSet->prepareNonLinear() is called
+
+ PinList m_pinList;
+ ElementList m_elementList;
+ ElementSet *m_elementSet;
+
+ //Stuff for caching
+ bool m_bCanCache;
+ LogicCacheNode * m_pLogicCacheBase;
+ unsigned m_logicOutCount;
+ LogicOut ** m_pLogicOut;
+
+ bool m_bCanAddChanged;
+ Circuit * m_pNextChanged[2];
+};
+
+#endif
diff --git a/src/electronics/simulation/currentsignal.cpp b/src/electronics/simulation/currentsignal.cpp
new file mode 100644
index 0000000..5e5388d
--- /dev/null
+++ b/src/electronics/simulation/currentsignal.cpp
@@ -0,0 +1,67 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "currentsignal.h"
+#include "element.h"
+#include "matrix.h"
+
+CurrentSignal::CurrentSignal( const double delta, const double current )
+ : Reactive::Reactive(delta)
+{
+ m_current = current;
+ m_oldCurrent = m_newCurrent = 0.;
+ m_numCNodes = 2;
+}
+
+CurrentSignal::~CurrentSignal()
+{
+}
+
+void CurrentSignal::setCurrent( const double i )
+{
+ if ( m_oldCurrent != m_newCurrent ) add_initial_dc();
+ m_newCurrent *= i/m_current; // Instead of calling step again, we can just "adjust" what the current should be
+ m_current = i;
+ add_initial_dc();
+}
+
+
+void CurrentSignal::add_map()
+{
+ // We don't need a map for current signal :-)
+}
+
+
+void CurrentSignal::add_initial_dc()
+{
+ if ( !b_status )
+ return;
+
+ if ( m_newCurrent == m_oldCurrent )
+ return;
+
+ b_i( 0 ) -= m_newCurrent-m_oldCurrent;
+ b_i( 1 ) += m_newCurrent-m_oldCurrent;
+
+ m_oldCurrent = m_newCurrent;
+}
+
+void CurrentSignal::updateCurrents()
+{
+ m_cnodeI[1] = m_newCurrent;
+ m_cnodeI[0] = -m_newCurrent;
+}
+
+void CurrentSignal::time_step()
+{
+ add_initial_dc(); // Make sure our old and new are synced
+ m_newCurrent = m_current*advance();
+ add_initial_dc();
+}
diff --git a/src/electronics/simulation/currentsignal.h b/src/electronics/simulation/currentsignal.h
new file mode 100644
index 0000000..b217b79
--- /dev/null
+++ b/src/electronics/simulation/currentsignal.h
@@ -0,0 +1,43 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef CURRENTSIGNAL_H
+#define CURRENTSIGNAL_H
+
+#include "reactive.h"
+#include "elementsignal.h"
+
+/**
+@short CurrentSignal
+@author David saxton
+*/
+class CurrentSignal : public Reactive, public ElementSignal
+{
+public:
+ CurrentSignal( const double delta, const double current );
+ virtual ~CurrentSignal();
+
+ virtual Element::Type type() const { return Element_CurrentSignal; }
+ void setCurrent( const double current );
+ double current() { return m_current; }
+ virtual void time_step();
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_current; // Current
+ double m_oldCurrent; // Old calculated current
+ double m_newCurrent; // New calculated current
+};
+
+#endif
diff --git a/src/electronics/simulation/currentsource.cpp b/src/electronics/simulation/currentsource.cpp
new file mode 100644
index 0000000..675b0b7
--- /dev/null
+++ b/src/electronics/simulation/currentsource.cpp
@@ -0,0 +1,64 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "currentsource.h"
+#include "elementset.h"
+#include "matrix.h"
+
+CurrentSource::CurrentSource( const double current )
+ : Element::Element()
+{
+ m_i = current;
+ m_numCNodes = 2;
+}
+
+
+CurrentSource::~CurrentSource()
+{
+}
+
+
+void CurrentSource::setCurrent( const double i )
+{
+ if ( i == m_i ) return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ // Remove the old current
+ m_i = -m_i;
+ add_initial_dc();
+
+ m_i = i;
+ add_initial_dc();
+}
+
+
+void CurrentSource::add_map()
+{
+ // We don't need a map for current source :-)
+}
+
+
+void CurrentSource::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ b_i( 0 ) -= m_i;
+ b_i( 1 ) += m_i;
+}
+
+
+void CurrentSource::updateCurrents()
+{
+ m_cnodeI[0] = -m_i;
+ m_cnodeI[1] = m_i;
+}
diff --git a/src/electronics/simulation/currentsource.h b/src/electronics/simulation/currentsource.h
new file mode 100644
index 0000000..9b72e0d
--- /dev/null
+++ b/src/electronics/simulation/currentsource.h
@@ -0,0 +1,39 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef CURRENTSOURCE_H
+#define CURRENTSOURCE_H
+
+#include "element.h"
+
+/**
+cnode n0 has current flowing otu of it, cnode n1 has current flowing into it
+@author David Saxton
+@short Current Source
+*/
+class CurrentSource : public Element
+{
+public:
+ CurrentSource( const double current );
+ virtual ~CurrentSource();
+
+ virtual Type type() const { return Element_CurrentSource; }
+ void setCurrent( const double i );
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_i; // Current
+};
+
+#endif
diff --git a/src/electronics/simulation/diode.cpp b/src/electronics/simulation/diode.cpp
new file mode 100644
index 0000000..e13d478
--- /dev/null
+++ b/src/electronics/simulation/diode.cpp
@@ -0,0 +1,198 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include <vector>
+#include "diode.h"
+#include "elementset.h"
+#include "matrix.h"
+
+#include <cmath>
+
+
+//BEGIN class Diode Settings
+DiodeSettings::DiodeSettings()
+{
+ reset();
+}
+
+
+void DiodeSettings::reset()
+{
+ I_S = 1e-15;
+ N = 1.0;
+ V_B = 4.7;
+// R = 0.001;
+}
+//END class Diode Settings
+
+
+
+//BEGIN class Diode
+Diode::Diode()
+ : NonLinear()
+{
+ m_numCNodes = 2;
+ g_new = g_old = I_new = I_old = V_prev = 0.0;
+}
+
+
+Diode::~Diode()
+{
+}
+
+
+void Diode::add_map()
+{
+ if (!b_status)
+ return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[0]->n(), Map::et_unstable, false );
+ }
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse( p_cnode[1]->n(), p_cnode[1]->n(), Map::et_unstable, false );
+ }
+
+ if ( !p_cnode[0]->isGround && !p_cnode[1]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[1]->n(), Map::et_unstable, false );
+ p_A->setUse( p_cnode[1]->n(), p_cnode[0]->n(), Map::et_unstable, false );
+ }
+}
+
+
+void Diode::add_initial_dc()
+{
+ g_new = g_old = I_new = I_old = V_prev = 0.0;
+ update_dc();
+}
+
+
+double Diode::current() const
+{
+ if (!b_status)
+ return 0.0;
+
+ double I;
+ calcIg( p_cnode[0]->v - p_cnode[1]->v, & I, 0 );
+
+ return I;
+}
+
+
+void Diode::updateCurrents()
+{
+ if (!b_status)
+ return;
+
+ m_cnodeI[1] = current();
+ m_cnodeI[0] = -m_cnodeI[1];
+}
+
+
+void Diode::update_dc()
+{
+ if (!b_status)
+ return;
+
+ calc_eq();
+
+ A_g( 0, 0 ) += g_new - g_old;
+ A_g( 1, 1 ) += g_new - g_old;
+ A_g( 0, 1 ) -= g_new - g_old;
+ A_g( 1, 0 ) -= g_new - g_old;
+
+ b_i( 0 ) -= I_new - I_old;
+ b_i( 1 ) += I_new - I_old;
+
+ g_old = g_new;
+ I_old = I_new;
+}
+
+
+
+#ifndef MIN
+# define MIN(x,y) (((x) < (y)) ? (x) : (y))
+#endif
+
+
+
+void Diode::calc_eq()
+{
+ double I_S = m_diodeSettings.I_S;
+ double N = m_diodeSettings.N;
+ double V_B = m_diodeSettings.V_B;
+// double R = m_diodeSettings.R;
+
+ double v = p_cnode[0]->v - p_cnode[1]->v;
+
+ // adjust voltage to help convergence
+ double V_crit = diodeCriticalVoltage( I_S, N * V_T );
+ if (V_B != 0 && v < MIN (0, -V_B + 10 * N * V_T))
+ {
+ double V = -(v + V_B);
+ V = diodeVoltage( V, -(V_prev + V_B), V_T * N, V_crit );
+ v = -(V + V_B);
+ }
+ else
+ v = diodeVoltage( v, V_prev, V_T * N, V_crit );
+
+ V_prev = v;
+
+ double I_D;
+ calcIg( v, & I_D, & g_new );
+
+ I_new = I_D - (v * g_new);
+}
+
+
+void Diode::calcIg( double V, double * I_D, double * g ) const
+{
+ double I_S = m_diodeSettings.I_S;
+ double N = m_diodeSettings.N;
+ double V_B = m_diodeSettings.V_B;
+// double R = m_diodeSettings.R;
+
+ double gtiny = (V < - 10 * V_T * N && V_B != 0) ? I_S : 0;
+
+ if ( V >= (-3 * N * V_T) )
+ {
+ if ( g )
+ *g = diodeConductance( V, I_S, V_T * N ) + gtiny;
+ *I_D = diodeCurrent( V, I_S, V_T * N ) + (gtiny * V);
+ }
+ else if ( V_B == 0 || V >= -V_B )
+ {
+ double a = (3 * N * V_T) / (V * M_E);
+ a = a * a * a;
+ *I_D = (-I_S * (1 + a)) + (gtiny * V);
+ if ( g )
+ *g = ((I_S * 3 * a) / V) + gtiny;
+ }
+ else
+ {
+ double a = exp( -(V_B + V) / N / V_T );
+ *I_D = (-I_S * a) + (gtiny * V);
+ if ( g )
+ *g = I_S * a / V_T / N + gtiny;
+ }
+}
+
+
+void Diode::setDiodeSettings( const DiodeSettings & settings )
+{
+ m_diodeSettings = settings;
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+}
+//END class Diode
+
diff --git a/src/electronics/simulation/diode.h b/src/electronics/simulation/diode.h
new file mode 100644
index 0000000..0b13946
--- /dev/null
+++ b/src/electronics/simulation/diode.h
@@ -0,0 +1,73 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef DIODE_H
+#define DIODE_H
+
+#include "nonlinear.h"
+
+class DiodeSettings
+{
+ public:
+ DiodeSettings();
+ void reset();
+
+ double I_S; ///< Diode saturation current
+ double N; ///< Emission coefficient
+ double V_B; ///< Reverse breakdown
+// double R; ///< Series resistance
+};
+
+
+/**
+This simulates a diode. The simulated diode characteristics are:
+
+@li I_s: Diode saturation current
+@li V_T: Thermal voltage = kT/4 = 25mV at 20 C
+@li n: Emission coefficient, typically between 1 and 2
+@li V_RB: Reverse breakdown (large negative voltage)
+@li G_RB: Reverse breakdown conductance
+@li R_D: Base resistance of diode
+
+@short Simulates the electrical property of diode-ness
+@author David Saxton
+*/
+class Diode : public NonLinear
+{
+ public:
+ Diode();
+ virtual ~Diode();
+
+ virtual void update_dc();
+ virtual void add_initial_dc();
+ virtual void add_map();
+ virtual Element::Type type() const { return Element_Diode; }
+ DiodeSettings settings() const { return m_diodeSettings; }
+ void setDiodeSettings( const DiodeSettings & settings );
+ /**
+ * Returns the current flowing through the diode
+ */
+ double current() const;
+
+ protected:
+ virtual void updateCurrents();
+ void calc_eq();
+
+ void calcIg( double V, double * I, double * g ) const;
+
+ double g_new, g_old;
+ double I_new, I_old;
+ double V_prev;
+
+ DiodeSettings m_diodeSettings;
+};
+
+#endif
+
diff --git a/src/electronics/simulation/element.cpp b/src/electronics/simulation/element.cpp
new file mode 100644
index 0000000..2411897
--- /dev/null
+++ b/src/electronics/simulation/element.cpp
@@ -0,0 +1,193 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "element.h"
+#include "elementset.h"
+
+#include <assert.h>
+
+#include <kdebug.h>
+
+Element::Element()
+{
+ b_status = false;
+ p_A = 0l;
+ p_eSet = 0l;
+ p_b = 0l;
+ b_componentDeleted = false;
+ b_eSetDeleted = true;
+
+ for ( int i=0; i<8; i++ )
+ p_cnode[i] = 0l;
+
+ resetCurrents();
+
+ for ( int i=0; i<4; i++ )
+ p_cbranch[i] = 0l;
+
+ m_numCBranches = 0;
+ m_numCNodes = 0;
+}
+
+Element::~ Element()
+{
+}
+
+void Element::resetCurrents()
+{
+ for ( int i=0; i<8; i++ )
+ m_cnodeI[i] = 0.0;
+}
+
+void Element::setElementSet( ElementSet *c )
+{
+ assert(!b_componentDeleted);
+ assert(!p_eSet);
+ if (!c) return elementSetDeleted();
+ b_eSetDeleted = false;
+ p_eSet = c;
+ p_A = p_eSet->matrix();
+ p_b = p_eSet->b();
+ updateStatus();
+}
+
+void Element::componentDeleted()
+{
+// assert(!b_componentDeleted);
+ if (b_componentDeleted)
+ {
+ // Something strange happened here....
+ }
+ if (b_eSetDeleted) return delete this;
+ b_componentDeleted = true;
+ b_status = false;
+// kdDebug() << "Element::componentDeleted(): Setting b_status to false, this="<<this<<endl;
+
+ p_eSet = 0l;
+ p_A = 0l;
+ p_b = 0l;
+ setCNodes();
+ setCBranches();
+}
+
+void Element::elementSetDeleted()
+{
+// assert(!b_eSetDeleted);
+ if (b_eSetDeleted)
+ {
+ // Something strange happened here....
+ }
+ if (b_componentDeleted) return delete this;
+ b_eSetDeleted = true;
+ b_status = false;
+// kdDebug() << "Element::elementSetDeleted(): Setting b_status to false, this="<<this<<endl;
+
+ p_eSet = 0l;
+ p_A = 0l;
+ p_b = 0l;
+ setCNodes();
+ setCBranches();
+}
+
+
+void Element::setCNodes( const int n0, const int n1, const int n2, const int n3 )
+{
+ if ( !p_eSet )
+ {
+// cerr << "Element::setCNodes: can't set nodes without circuit!"<<endl;
+ for ( int i=0; i<8; i++ )
+ p_cnode[i] = 0l;
+ return;
+ }
+
+ // MAX_CNODES-1 should match the last array index below.
+ assert( MAX_CNODES == 4 );
+ p_cnode[0] = (n0>-1)?p_eSet->cnodes()[n0]:(n0==-1?p_eSet->ground():0l);
+ p_cnode[1] = (n1>-1)?p_eSet->cnodes()[n1]:(n1==-1?p_eSet->ground():0l);
+ p_cnode[2] = (n2>-1)?p_eSet->cnodes()[n2]:(n2==-1?p_eSet->ground():0l);
+ p_cnode[3] = (n3>-1)?p_eSet->cnodes()[n3]:(n3==-1?p_eSet->ground():0l);
+ updateStatus();
+}
+
+void Element::setCBranches( const int b0, const int b1, const int b2, const int b3 )
+{
+ if ( !p_eSet )
+ {
+// cerr << "Element::setCBranches: can't set branches without circuit!"<<endl;
+ for ( int i=0; i<4; i++ ) p_cbranch[i] = 0l;
+ return;
+ }
+ p_cbranch[0] = (b0>-1)?p_eSet->cbranches()[b0]:0l;
+ p_cbranch[1] = (b1>-1)?p_eSet->cbranches()[b1]:0l;
+ p_cbranch[2] = (b2>-1)?p_eSet->cbranches()[b2]:0l;
+ p_cbranch[3] = (b3>-1)?p_eSet->cbranches()[b3]:0l;
+ updateStatus();
+}
+
+bool Element::updateStatus()
+{
+ // First, set status to false if all nodes in use are ground
+ b_status = false;
+ for ( int i=0; i<m_numCNodes; i++ )
+ {
+ b_status |= p_cnode[i]?!p_cnode[i]->isGround:false;
+ }
+
+ // Set status to false if any of the nodes are not set
+ for ( int i=0; i<m_numCNodes; i++ )
+ {
+ if (!p_cnode[i]) b_status = false;
+ }
+
+ // Finally, set status to false if not all the required branches are set
+ for ( int i=0; i<m_numCBranches; i++ )
+ {
+ if (!p_cbranch[i]) b_status = false;
+ }
+
+ // Finally, check for various pointers
+ if ( !p_eSet || !p_A || !p_b ) b_status = false;
+
+ if (!b_status)
+ {
+ resetCurrents();
+ }
+ // And return the status :-)
+// kdDebug() << "Element::updateStatus(): Setting b_status to "<<(b_status?"true":"false")<<" this="<<this<<endl;
+ return b_status;
+}
+
+double Element::cbranchCurrent( const int branch )
+{
+ if ( !b_status || branch<0 || branch>=m_numCBranches ) return 0.;
+ return (*p_cbranch)[branch].i;
+}
+
+double Element::cnodeVoltage( const int node )
+{
+ if ( !b_status || node<0 || node>=m_numCNodes ) return 0.;
+ return (*p_cnode)[node].v;
+}
+
+
+CNode::CNode()
+{
+ m_n = 0;
+ v = 0.0;
+ isGround = false;
+}
+
+CBranch::CBranch()
+{
+ m_n = 0;
+ i = 0.0;
+}
+
+
diff --git a/src/electronics/simulation/element.h b/src/electronics/simulation/element.h
new file mode 100644
index 0000000..e05de46
--- /dev/null
+++ b/src/electronics/simulation/element.h
@@ -0,0 +1,255 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef ELEMENT_H
+#define ELEMENT_H
+
+#include "elementset.h"
+#include "matrix.h"
+
+class ElementSet;
+class Vector;
+typedef unsigned int uint;
+
+const double T = 300.; // Temperature in Kelvin
+const double K = 1.38e-23; // Boltzmann's constant
+const double q = 1.602e-19; // Charge on an electron
+const double V_T = K*T/q; // Thermal voltage
+const double gmin = 1e-12; // Minimum parallel conductance used in dc domain
+
+class CNode
+{
+public:
+ CNode();
+ void set_n( const uint n ) { m_n=n; }
+ uint n() { return m_n; }
+ double v; // Voltage on node. This is set from the last calculated voltage.
+ bool isGround; // True for ground nodes. Obviously, you should ignore n and v if this is true
+private:
+ uint m_n; // CNode number
+};
+
+class CBranch
+{
+public:
+ CBranch();
+ void set_n( const uint n ) { m_n=n; }
+ uint n() { return m_n; }
+ double i; // Current flowing through branch. This is set from the last calculated current.
+private:
+ uint m_n; // CBranch number
+};
+
+const int MAX_CNODES = 4;
+
+// Default node number that represents no node (remember that
+// Ground node is -1, and the rest are numbered from 0 to n-1
+const int noCNode = -2;
+// Likewise for branch (although there is no "ground" branch;
+// it is merely -2 for likeness with noCNode)
+const int noBranch = -2;
+
+/**
+@short Represents a circuit element (such as resistance)
+@author David Saxton
+*/
+class Element
+{
+public:
+ enum Type
+ {
+ Element_BJT,
+ Element_Capacitance,
+ Element_CCCS,
+ Element_CCVS,
+ Element_CurrentSignal,
+ Element_CurrentSource,
+ Element_Diode,
+ Element_Inductance,
+ Element_LogicIn,
+ Element_LogicOut,
+ Element_OpAmp,
+ Element_Resistance,
+ Element_VCCS,
+ Element_VCVS,
+ Element_VoltagePoint,
+ Element_VoltageSignal,
+ Element_VoltageSource
+ };
+
+ Element();
+ virtual ~Element();
+ /**
+ * This must be called when the circuit is changed. The function will get
+ * all the required pointers from ElementSet
+ */
+ virtual void setElementSet( ElementSet *c );
+ /**
+ * Returns a pointer to the current element set
+ */
+ ElementSet *elementSet() { return p_eSet; }
+ /**
+ * Tells the element which nodes to use. Remember that -1 is ground. You
+ * should refer to the individual elements for which nodes are used for what.
+ */
+ void setCNodes( const int n0 = noCNode, const int n1 = noCNode, const int n2 = noCNode, const int n3 = noCNode );
+ /**
+ * Tells the element it's branch numbers (if it should have one). Not
+ * all elements use this.
+ */
+ void setCBranches( const int b0 = noBranch, const int b1 = noBranch, const int b2 = noBranch, const int b3 = noBranch );
+ /**
+ * Returns a pointer to the given CNode
+ */
+ CNode *cnode( const uint num ) { return p_cnode[num]; }
+ /**
+ * Returns a pointer to the given CNode
+ */
+ CBranch *cbranch( const uint num ) { return p_cbranch[num]; }
+ /**
+ * Returns the number of branches used by the element
+ */
+ int numCBranches() { return m_numCBranches; }
+ /**
+ * Returns the number of circuit nodes used by the element
+ */
+ int numCNodes() { return m_numCNodes; }
+ /**
+ * Call this function to tell the element to calculate the
+ * current flowing *into* it's cnodes *from* the element. You
+ * can get the currents with m_cnodeI. Child class must implement this function.
+ */
+ virtual void updateCurrents() = 0;
+ /**
+ * Returns true for reactive elements that need stepping for numerical-integration
+ * (such as capacitors)
+ */
+ virtual bool isReactive() { return false; }
+ /**
+ * Returns true for NonLinear elements that need iteration to converge to a solution
+ * as the matrix A is a function of x.
+ */
+ virtual bool isNonLinear() { return false; }
+ /**
+ * Returns the type of element
+ */
+ virtual Type type() const = 0;
+ /**
+ * Call this function to tell the element to add its map to the matrix in use
+ */
+ virtual void add_map() {};
+ /**
+ * Does the required MNA stuff. This should be called from ElementSet when necessary.
+ */
+ virtual void add_initial_dc() = 0;
+ /**
+ * This is called from the Component destructor. When elementSetDeleted has
+ * also been called, this class will delete itself.
+ */
+ void componentDeleted();
+ void elementSetDeleted();
+
+ double m_cnodeI[8]; ///< Current flowing into the cnodes from the element
+ double cbranchCurrent( const int branch );
+ double cnodeVoltage( const int node );
+
+protected:
+ /**
+ * Resets all calculated currents in the nodes to 0
+ */
+ void resetCurrents();
+
+ inline double & A_g( uint i, uint j );
+ inline double & A_b( uint i, uint j );
+ inline double & A_c( uint i, uint j );
+ inline double & A_d( uint i, uint j );
+
+ inline double & b_i( uint i );
+ inline double & b_v( uint i );
+
+ ElementSet *p_eSet;
+ Matrix *p_A;
+ Vector *p_b;
+ CNode *p_cnode[MAX_CNODES];
+ CBranch *p_cbranch[4];
+
+ /**
+ * True when the element can do add_initial_dc(), i.e. when it has
+ * pointers to the circuit, and at least one of its nodes is not ground.
+ */
+ bool b_status;
+ /**
+ * Update the status, returning b_status
+ */
+ virtual bool updateStatus();
+ /**
+ * Set by child class - the number of branches that the element uses
+ * Typically, this is 0, but could be 1 (e.g. independent voltage source)
+ * or 2 (e.g. cccs)
+ */
+ int m_numCBranches;
+ /**
+ * Set by child class - the number of circuit nodes that the element uses
+ */
+ int m_numCNodes;
+
+private:
+ bool b_componentDeleted;
+ bool b_eSetDeleted;
+ double m_temp;
+};
+
+
+double & Element::A_g( uint i, uint j )
+{
+ if ( p_cnode[i]->isGround || p_cnode[j]->isGround )
+ return m_temp;
+ return p_A->g( p_cnode[i]->n(), p_cnode[j]->n() );
+}
+
+
+double & Element::A_b( uint i, uint j )
+{
+ if ( p_cnode[i]->isGround )
+ return m_temp;
+ return p_A->b( p_cnode[i]->n(), p_cbranch[j]->n() );
+}
+
+
+double & Element::A_c( uint i, uint j )
+{
+ if ( p_cnode[j]->isGround )
+ return m_temp;
+ return p_A->c( p_cbranch[i]->n(), p_cnode[j]->n() );
+}
+
+
+double & Element::A_d( uint i, uint j )
+{
+ return p_A->d( p_cbranch[i]->n(), p_cbranch[j]->n() );
+}
+
+
+
+double & Element::b_i( uint i )
+{
+ if ( p_cnode[i]->isGround )
+ return m_temp;
+
+ return (*p_b)[ p_cnode[i]->n() ];
+}
+
+
+double & Element::b_v( uint i )
+{
+ return (*p_b)[ p_eSet->cnodeCount() + p_cbranch[i]->n() ];
+}
+
+#endif
diff --git a/src/electronics/simulation/elementset.cpp b/src/electronics/simulation/elementset.cpp
new file mode 100644
index 0000000..25057c2
--- /dev/null
+++ b/src/electronics/simulation/elementset.cpp
@@ -0,0 +1,253 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "bjt.h"
+#include "circuit.h"
+#include "elementset.h"
+#include "element.h"
+#include "logic.h"
+#include "matrix.h"
+#include "nonlinear.h"
+#include "vec.h"
+
+#include <kdebug.h>
+
+#include <cmath>
+#include <iostream>
+#include <assert.h>
+
+ElementSet::ElementSet( Circuit * circuit, const int n, const int m )
+{
+ m_pCircuit = circuit;
+ m_cn = n;
+ m_cb = m;
+ p_logicIn = 0l;
+ p_A = new Matrix( m_cn, m_cb );
+ p_b = new Vector(m_cn+m_cb);
+ p_x = new Vector(m_cn+m_cb);
+ p_x_prev = new Vector(m_cn+m_cb);
+ m_cbranches = new CBranch*[m_cb];
+ m_cnodes = new CNode*[m_cn];
+ for ( uint i=0; i<m_cn; i++ )
+ {
+ m_cnodes[i] = new CNode();
+ m_cnodes[i]->set_n(i);
+ }
+ for ( uint i=0; i<m_cb; i++ )
+ {
+ m_cbranches[i] = new CBranch();
+ m_cbranches[i]->set_n(i);
+ }
+ m_ground = new CNode();
+ m_ground->isGround = true;
+ b_containsNonLinear = false;
+}
+
+ElementSet::~ElementSet()
+{
+ const ElementList::iterator end = m_elementList.end();
+ for ( ElementList::iterator it = m_elementList.begin(); it != end; ++it )
+ {
+ // Note: By calling setElementSet(0l), we might have deleted it (the Element will commit
+ // suicide when both the the ElementSet and Component to which it belongs have deleted
+ // themselves). So be very careful it you plan to do anything with the (*it) pointer
+ if (*it) (*it)->elementSetDeleted();
+ }
+ for ( uint i=0; i<m_cn; i++ )
+ {
+ delete m_cnodes[i];
+ }
+ for ( uint i=0; i<m_cb; i++ )
+ {
+ delete m_cbranches[i];
+ }
+ delete[] m_cbranches;
+ delete[] m_cnodes;
+ delete[] p_logicIn;
+ delete m_ground;
+ delete p_A;
+ delete p_b;
+ delete p_x;
+ delete p_x_prev;
+}
+
+
+void ElementSet::setCacheInvalidated()
+{
+ m_pCircuit->setCacheInvalidated();
+}
+
+
+void ElementSet::addElement( Element *e )
+{
+ if ( !e || m_elementList.contains(e) ) return;
+ e->setElementSet(this);
+ m_elementList.append(e);
+ if ( e->isNonLinear() )
+ {
+ b_containsNonLinear = true;
+ m_cnonLinearList.append( static_cast<NonLinear*>(e) );
+ }
+}
+
+
+void ElementSet::createMatrixMap()
+{
+ p_A->createMap();
+
+
+ // And do our logic as well...
+
+ m_clogic = 0;
+ ElementList::iterator end = m_elementList.end();
+ for ( ElementList::iterator it = m_elementList.begin(); it != end; ++it )
+ {
+ if ( dynamic_cast<LogicIn*>(*it) )
+ m_clogic++;
+ }
+
+ p_logicIn = new LogicIn*[m_clogic];
+ int i=0;
+ for ( ElementList::iterator it = m_elementList.begin(); it != end; ++it )
+ {
+ if ( LogicIn * in = dynamic_cast<LogicIn*>(*it) )
+ p_logicIn[i++] = in;
+ }
+}
+
+
+void ElementSet::doNonLinear( int maxIterations, double maxErrorV, double maxErrorI )
+{
+// p_x_prev->reset();
+
+ // And now tell the cnodes and cbranches about their new voltages & currents
+ updateInfo();
+
+ const NonLinearList::iterator end = m_cnonLinearList.end();
+
+ int k = 0;
+ do
+ {
+ // Tell the nonlinear elements to update its J, A and b from the newly calculated x
+ for ( NonLinearList::iterator it = m_cnonLinearList.begin(); it != end; ++it )
+ (*it)->update_dc();
+
+ *p_x = *p_b;
+ p_A->performLU();
+ p_A->fbSub(p_x);
+ updateInfo();
+
+ // Now, check for convergence
+ bool converged = true;
+ for ( unsigned i = 0; i < m_cn; ++i )
+ {
+ double diff = std::abs( (*p_x_prev)[i] - (*p_x)[i] );
+ if ( diff > maxErrorI )
+ {
+ converged = false;
+ break;
+ }
+ }
+ if ( converged )
+ {
+ for ( unsigned i = m_cn; i < m_cn+m_cb; ++i )
+ {
+ double diff = std::abs( (*p_x_prev)[i] - (*p_x)[i] );
+ if ( diff > maxErrorV )
+ {
+ converged = false;
+ break;
+ }
+ }
+ }
+
+ *p_x_prev = *p_x;
+
+ if ( converged )
+ break;
+ }
+ while ( ++k < maxIterations );
+}
+
+
+bool ElementSet::doLinear( bool performLU )
+{
+ if ( b_containsNonLinear || (!p_b->isChanged() && ((performLU && !p_A->isChanged()) || !performLU)) )
+ return false;
+
+ if (performLU)
+ p_A->performLU();
+
+ *p_x = *p_b;
+ p_A->fbSub(p_x);
+ updateInfo();
+ p_b->setUnchanged();
+
+ return true;
+}
+
+
+void ElementSet::updateInfo()
+{
+ for ( uint i=0; i<m_cn; i++ )
+ {
+ const double v = (*p_x)[i];
+ if (std::isfinite(v)) {
+ m_cnodes[i]->v = v;
+ } else {
+ (*p_x)[i] = 0.;
+ m_cnodes[i]->v = 0.;
+ }
+ }
+ for ( uint i=0; i<m_cb; i++ )
+ {
+ // NOTE: I've used lowercase and uppercase "I" here, so be careful!
+ const double I = (*p_x)[i+m_cn];
+ if (std::isfinite(I)) {
+ m_cbranches[i]->i = I;
+ } else {
+ (*p_x)[i+m_cn] = 0.;
+ m_cbranches[i]->i = 0.;
+ }
+ }
+
+ // Tell logic to check themselves
+ for ( uint i=0; i<m_clogic; ++i )
+ {
+ p_logicIn[i]->check();
+ }
+}
+
+
+void ElementSet::displayEquations()
+{
+ std::cout.setf(std::ios_base::fixed);
+ std::cout.precision(5);
+ std::cout.setf(std::ios_base::showpoint);
+ std::cout << "A x = b :"<<std::endl;
+ for ( uint i=0; i<m_cn+m_cb; i++ )
+ {
+ std::cout << "( ";
+ for ( uint j=0; j<m_cn+m_cb; j++ )
+ {
+ const double value = p_A->g(i,j);
+// if ( value > 0 ) cout <<"+";
+// else if ( value == 0 ) cout <<" ";
+ std::cout.width(10);
+ std::cout << value<<" ";
+ }
+ std::cout << ") ( "<<(*p_x)[i]<<" ) = ( ";
+ std::cout<<(*p_b)[i]<<" )"<<std::endl;
+ }
+ std::cout << "A_LU:"<<std::endl;
+ p_A->displayLU();
+}
+
+
diff --git a/src/electronics/simulation/elementset.h b/src/electronics/simulation/elementset.h
new file mode 100644
index 0000000..c7ef7ca
--- /dev/null
+++ b/src/electronics/simulation/elementset.h
@@ -0,0 +1,131 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef ELEMENTSET_H
+#define ELEMENTSET_H
+
+#include <vector>
+
+#include <qvaluelist.h>
+
+class CBranch;
+class Circuit;
+class CNode;
+class Element;
+class ElementSet;
+class LogicIn;
+class Matrix;
+class NonLinear;
+class Vector;
+
+typedef QValueList<Element*> ElementList;
+typedef QValueList<NonLinear*> NonLinearList;
+
+/**
+Steps in simulation of a set of elements:
+(1) Create this class with given number of nodes "n" and voltage sources "m"
+(2) Add the various elements with addElement.
+(3) Call performDC()
+(4) Get the nodal voltages and voltage currents with x()
+(5) Repeat steps 3 and 4 if necessary for further transient analysis.
+
+This class shouldn't be confused with the Circuit class, but considered a helper class to Circuit.
+Circuit will handle the simulation of a set of components over time. This just finds the DC-operating
+point of the circuit for a given set of elements.
+
+@short Handles a set of circuit elements
+@author David Saxton
+*/
+class ElementSet
+{
+public:
+ /**
+ * Create a new circuit, with "n" nodes and "m" voltage sources.
+ * After creating the circuit, you must call setGround to specify
+ * the ground nodes, before adding any elements.
+ */
+ ElementSet( Circuit * circuit, const int n, const int m );
+ /**
+ * Destructor. Note that only the matrix and supporting data is deleted.
+ * i.e. Any elements added to the circuit will not be deleted.
+ */
+ ~ElementSet();
+ Circuit * circuit() const { return m_pCircuit; }
+ void addElement( Element *e );
+ void setCacheInvalidated();
+ /**
+ * Returns the matrix in use. This is created once on the creation of the ElementSet
+ * class, and deleted in the destructor, so the pointer returned will never change.
+ */
+ Matrix *matrix() const { return p_A; }
+ /**
+ * Returns the vector for b (i.e. the independent currents & voltages)
+ */
+ Vector *b() const { return p_b; }
+ /**
+ * Returns the vector for x (i.e. the currents & voltages at the branches and nodes)
+ */
+ Vector *x() const { return p_x; }
+ /**
+ * @return if we have any nonlinear elements (e.g. diodes, tranaistors).
+ */
+ bool containsNonLinear() const { return b_containsNonLinear; }
+ /**
+ * Solves for nonlinear elements, or just does linear if it doesn't contain
+ * any nonlinear.
+ */
+ void doNonLinear( int maxIterations, double maxErrorV = 1e-9, double maxErrorI = 1e-12 );
+ /**
+ * Solves for linear and logic elements.
+ * @returns true if anything changed
+ */
+ bool doLinear( bool performLU );
+ CBranch **cbranches() const { return m_cbranches; }
+ CNode **cnodes() const { return m_cnodes; }
+ CNode *ground() const { return m_ground; }
+ /**
+ * Returns the number of nodes in the circuit (excluding ground 'nodes')
+ */
+ int cnodeCount() const { return m_cn; }
+ /**
+ * Returns the number of voltage sources in the circuit
+ */
+ int cbranchCount() const { return m_cb; }
+
+ void createMatrixMap();
+ /**
+ * Displays the matrix equations Ax=b and J(dx)=-r
+ */
+ void displayEquations();
+ /**
+ * Update the nodal voltages and branch currents from the x vector
+ */
+ void updateInfo();
+
+private:
+ Matrix *p_A;
+ Vector *p_x;
+ Vector *p_x_prev;
+ Vector *p_b;
+ uint m_cn;
+ uint m_cb;
+ ElementList m_elementList;
+ NonLinearList m_cnonLinearList;
+ CBranch **m_cbranches; // Pointer to an array of cbranches
+ CNode **m_cnodes; // Pointer to an array of cnodes
+ CNode *m_ground;
+ uint m_clogic;
+ LogicIn **p_logicIn;
+ bool b_containsNonLinear;
+ Circuit * m_pCircuit;
+};
+
+#endif
+
diff --git a/src/electronics/simulation/elementsignal.cpp b/src/electronics/simulation/elementsignal.cpp
new file mode 100644
index 0000000..31d7d78
--- /dev/null
+++ b/src/electronics/simulation/elementsignal.cpp
@@ -0,0 +1,64 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementsignal.h"
+#include <cmath>
+
+ElementSignal::ElementSignal()
+{
+ m_type = ElementSignal::st_sinusoidal;
+ m_time = 0.;
+ m_frequency = 0.;
+}
+
+ElementSignal::~ElementSignal()
+{
+}
+
+void ElementSignal::setStep( double delta, Type type, double frequency )
+{
+ m_type = type;
+ m_delta = delta;
+ m_frequency = frequency;
+ m_omega = 6.283185307179586*m_frequency;
+ m_time = 1./(4.*m_frequency);
+}
+
+double ElementSignal::advance()
+{
+ m_time += m_delta;
+ if ( m_time >= 1./m_frequency ) m_time -= 1./m_frequency;
+
+ switch (m_type)
+ {
+ case ElementSignal::st_sawtooth:
+ {
+ // TODO Sawtooth signal
+ return 0.;
+ }
+ case ElementSignal::st_square:
+ {
+ return (sin(m_time*m_omega)>=0)?1:-1;
+ }
+ case ElementSignal::st_triangular:
+ {
+ // TODO Triangular signal
+ return 0.;
+ }
+ case ElementSignal::st_sinusoidal:
+ default:
+ {
+ return sin(m_time*m_omega);
+ }
+ }
+}
+
+
+
diff --git a/src/electronics/simulation/elementsignal.h b/src/electronics/simulation/elementsignal.h
new file mode 100644
index 0000000..c26bea1
--- /dev/null
+++ b/src/electronics/simulation/elementsignal.h
@@ -0,0 +1,45 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef ELEMENTSIGNAL_H
+#define ELEMENTSIGNAL_H
+
+/**
+@short Provides different signals
+@author David Saxton
+*/
+class ElementSignal
+{
+public:
+ enum Type
+ {
+ st_sinusoidal,
+ st_square,
+ st_sawtooth,
+ st_triangular
+ };
+ ElementSignal();
+ ~ElementSignal();
+
+ void setStep( double delta, Type type, double frequency );
+ /**
+ * Advances the timer, returns amplitude (between -1 and 1)
+ */
+ double advance();
+
+protected:
+ Type m_type;
+ double m_time;
+ double m_frequency;
+ double m_delta;
+ double m_omega; // Used for sinusoidal signal
+};
+
+#endif
diff --git a/src/electronics/simulation/inductance.cpp b/src/electronics/simulation/inductance.cpp
new file mode 100644
index 0000000..22c5f9e
--- /dev/null
+++ b/src/electronics/simulation/inductance.cpp
@@ -0,0 +1,126 @@
+/***************************************************************************
+ * Copyright (C) 2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementset.h"
+#include "inductance.h"
+#include "matrix.h"
+
+Inductance::Inductance( double inductance, double delta )
+ : Reactive(delta)
+{
+ m_inductance = inductance;
+ r_eq_old = v_eq_old = 0.0;
+ m_numCNodes = 2;
+ m_numCBranches = 1;
+ setMethod( Inductance::m_euler );
+}
+
+
+Inductance::~Inductance()
+{
+}
+
+
+void Inductance::setInductance( double i )
+{
+ m_inductance = i;
+}
+
+
+void Inductance::add_initial_dc()
+{
+ A_c( 0, 0 ) = 1;
+ A_b( 0, 0 ) = 1;
+ A_c( 0, 1 ) = -1;
+ A_b( 1, 0 ) = -1;
+
+ // The adding of r_eg and v_eq will be done for us by time_step.
+ // So for now, just reset the constants used.
+ r_eq_old = v_eq_old = 0.0;
+}
+
+
+void Inductance::updateCurrents()
+{
+ if (!b_status)
+ return;
+ m_cnodeI[0] = p_cbranch[0]->i;
+ m_cnodeI[1] = -m_cnodeI[0];
+}
+
+
+void Inductance::add_map()
+{
+ if (!b_status)
+ return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_constant, true );
+ p_A->setUse_b( p_cnode[0]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ }
+
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[1]->n(), Map::et_constant, true );
+ p_A->setUse_b( p_cnode[1]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ }
+
+ p_A->setUse_d( p_cbranch[0]->n(), p_cbranch[0]->n(), Map::et_unstable, false );
+}
+
+
+void Inductance::time_step()
+{
+ if (!b_status) return;
+
+ double i = p_cbranch[0]->i;
+ double v_eq_new = 0.0, r_eq_new = 0.0;
+
+ if ( m_method == Inductance::m_euler )
+ {
+ r_eq_new = m_inductance/m_delta;
+ v_eq_new = -i*r_eq_new;
+ }
+ else if ( m_method == Inductance::m_trap ) {
+ // TODO Implement + test trapezoidal method
+ r_eq_new = 2.0*m_inductance/m_delta;
+ }
+
+ if ( r_eq_old != r_eq_new )
+ {
+ A_d( 0, 0 ) -= r_eq_new - r_eq_old;
+ }
+
+ if ( v_eq_new != v_eq_old )
+ {
+ b_v( 0 ) += v_eq_new - v_eq_old;
+ }
+
+ r_eq_old = r_eq_new;
+ v_eq_old = v_eq_new;
+}
+
+
+bool Inductance::updateStatus()
+{
+ b_status = Reactive::updateStatus();
+ if ( m_method == Inductance::m_none )
+ b_status = false;
+ return b_status;
+}
+
+
+void Inductance::setMethod( Method m )
+{
+ m_method = m;
+ updateStatus();
+}
+
diff --git a/src/electronics/simulation/inductance.h b/src/electronics/simulation/inductance.h
new file mode 100644
index 0000000..46ccb09
--- /dev/null
+++ b/src/electronics/simulation/inductance.h
@@ -0,0 +1,55 @@
+/***************************************************************************
+ * Copyright (C) 2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef INDUCTANCE_H
+#define INDUCTANCE_H
+
+#include "reactive.h"
+
+/**
+
+@author David Saxton
+*/
+class Inductance : public Reactive
+{
+ public:
+ enum Method
+ {
+ m_none, // None
+ m_euler, // Backward Euler
+ m_trap // Trapezoidal (currently unimplemented)
+ };
+ Inductance( double capacitance, double delta );
+ virtual ~Inductance();
+
+ virtual Type type() const { return Element_Inductance; }
+ /**
+ * Set the stepping use for numerical integration of inductance, and the
+ * interval between successive updates.
+ */
+ void setMethod( Method m );
+ virtual void time_step();
+ virtual void add_initial_dc();
+ void setInductance( double i );
+ virtual void add_map();
+
+ protected:
+ virtual void updateCurrents();
+ virtual bool updateStatus();
+
+ private:
+ double m_inductance; // Inductance
+ Method m_method; // Method of integration
+
+ double r_eq_old;
+ double v_eq_old;
+};
+
+#endif
diff --git a/src/electronics/simulation/logic.cpp b/src/electronics/simulation/logic.cpp
new file mode 100644
index 0000000..031dd2e
--- /dev/null
+++ b/src/electronics/simulation/logic.cpp
@@ -0,0 +1,319 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include <vector>
+#include "circuit.h"
+#include "elementset.h"
+#include "logic.h"
+#include "matrix.h"
+#include "simulator.h"
+#include "src/core/ktlconfig.h"
+
+LogicIn::LogicIn( LogicConfig config )
+ : Element::Element()
+{
+ m_config = config;
+ m_pCallbackFunction = 0l;
+ m_numCNodes = 1;
+ m_bLastState = false;
+ m_pNextLogic = 0l;
+ setLogic(getConfig());
+}
+
+LogicIn::~LogicIn()
+{
+ Simulator::self()->removeLogicInReferences(this);
+}
+
+
+void LogicIn::setCallback( CallbackClass * object, CallbackPtr func )
+{
+ m_pCallbackFunction = func;
+ m_pCallbackObject = object;
+}
+
+
+void LogicIn::check()
+{
+ if (!b_status)
+ return;
+
+ bool newState;
+ if (m_bLastState)
+ {
+ // Was high, will still be high unless voltage is less than falling trigger
+ newState = p_cnode[0]->v > m_config.fallingTrigger;
+ }
+ else
+ {
+ // Was low, will still be low unless voltage is more than rising trigger
+ newState = p_cnode[0]->v > m_config.risingTrigger;
+ }
+
+ if ( m_pCallbackFunction && (newState != m_bLastState) )
+ {
+ m_bLastState = newState;
+ (m_pCallbackObject->*m_pCallbackFunction)(newState);
+ }
+ m_bLastState = newState;
+}
+
+
+void LogicIn::setLogic( LogicConfig config )
+{
+ m_config = config;
+ check();
+}
+
+
+void LogicIn::setElementSet( ElementSet *c )
+{
+ if (c)
+ m_pNextLogic = 0l;
+ else
+ m_cnodeI[0] = 0.;
+
+ Element::setElementSet(c);
+}
+
+
+void LogicIn::add_initial_dc()
+{
+}
+
+
+void LogicIn::updateCurrents()
+{
+}
+
+
+LogicConfig LogicIn::getConfig()
+{
+ LogicConfig c;
+ c.risingTrigger = KTLConfig::logicRisingTrigger();
+ c.fallingTrigger = KTLConfig::logicFallingTrigger();
+ c.output = KTLConfig::logicOutputHigh();
+ c.highImpedance = KTLConfig::logicOutputHighImpedance();
+ c.lowImpedance = KTLConfig::logicOutputLowImpedance();
+ return c;
+}
+
+
+LogicOut::LogicOut( LogicConfig config, bool _high )
+ : LogicIn(config)
+{
+ m_bCanAddChanged = true;
+ m_bOutputHighConductanceConst = false;
+ m_bOutputLowConductanceConst = false;
+ m_bOutputHighVoltageConst = false;
+ m_pNextChanged[0] = m_pNextChanged[1] = 0l;
+ m_pSimulator = 0l;
+ m_bUseLogicChain = false;
+ b_state = false;
+ m_numCNodes = 1;
+ m_vHigh = m_gHigh = m_gLow = 0.0;
+ m_old_g_out = m_g_out = 0.0;
+ m_old_v_out = m_v_out = 0.0;
+ setHigh(_high);
+
+ // Although we already call this function in LogicIn's constructor, our
+ // virtual function will not have got called, so we have to call it again.
+ setLogic(getConfig());
+}
+
+LogicOut::~LogicOut()
+{
+ if (!m_pSimulator)
+ m_pSimulator = Simulator::self();
+
+ // Note that although this function will get called in the destructor of
+ // LogicIn, we must call it here as well as it needs to be called before
+ // removeLogicOutReferences(this) is called.
+ m_pSimulator->removeLogicInReferences(this);
+
+ m_pSimulator->removeLogicOutReferences(this);
+}
+
+
+void LogicOut::setUseLogicChain( bool use )
+{
+ if (!m_pSimulator)
+ m_pSimulator = Simulator::self();
+
+ m_bUseLogicChain = use;
+ if (use)
+ setElementSet(0l);
+}
+
+
+void LogicOut::setElementSet( ElementSet *c )
+{
+ if (!m_pSimulator)
+ m_pSimulator = Simulator::self();
+
+ if (c)
+ {
+ m_bUseLogicChain = false;
+ m_pNextChanged[0] = m_pNextChanged[1] = 0l;
+ }
+
+ // NOTE Make sure that the next two lines are the same as those in setHigh and setLogic
+ m_g_out = b_state ? m_gHigh : m_gLow;
+ m_v_out = b_state ? m_vHigh : 0.0;
+
+ LogicIn::setElementSet(c);
+}
+
+
+void LogicOut::setOutputHighConductance( double g )
+{
+ m_bOutputHighConductanceConst = true;
+ if ( g == m_gHigh )
+ return;
+ m_gHigh = g;
+ configChanged();
+}
+
+
+void LogicOut::setOutputLowConductance( double g )
+{
+ m_bOutputLowConductanceConst = true;
+ if ( g == m_gLow )
+ return;
+ m_gLow = g;
+ configChanged();
+}
+
+
+void LogicOut::setOutputHighVoltage( double v )
+{
+ m_bOutputHighVoltageConst = true;
+ if ( v == m_vHigh )
+ return;
+ m_vHigh = v;
+ configChanged();
+}
+
+
+void LogicOut::setLogic( LogicConfig config )
+{
+ m_config = config;
+
+ if (!m_bOutputHighConductanceConst)
+ m_gHigh = 1.0/config.highImpedance;
+
+ if (!m_bOutputLowConductanceConst)
+ m_gLow = (config.lowImpedance == 0.0) ? 0.0 : 1.0/config.lowImpedance;
+
+ if (!m_bOutputHighVoltageConst)
+ m_vHigh = config.output;
+
+ configChanged();
+}
+
+
+void LogicOut::configChanged()
+{
+ if (m_bUseLogicChain)
+ return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ // Re-add the DC stuff using the new values
+
+ m_old_g_out = m_g_out;
+ m_old_v_out = m_v_out;
+
+ // NOTE Make sure that the next two lines are the same as those in setElementSet and setHigh
+ m_g_out = b_state ? m_gHigh : m_gLow;
+ m_v_out = b_state ? m_vHigh : 0.0;
+
+ add_initial_dc();
+
+ m_old_g_out = 0.;
+ m_old_v_out = 0.;
+
+ check();
+}
+
+
+void LogicOut::add_map()
+{
+ if (!b_status) return;
+
+ p_A->setUse( p_cnode[0]->n(), p_cnode[0]->n(), Map::et_variable, false );
+}
+
+
+void LogicOut::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ A_g( 0, 0 ) += m_g_out-m_old_g_out;
+ b_i( 0 ) += m_g_out*m_v_out-m_old_g_out*m_old_v_out;
+}
+
+void LogicOut::updateCurrents()
+{
+ if (m_bUseLogicChain)
+ {
+ m_cnodeI[0] = 0.;
+ return;
+ }
+ if (!b_status) {
+ return;
+ }
+ m_cnodeI[0] = (p_cnode[0]->v-m_v_out)*m_g_out;
+}
+
+void LogicOut::setHigh( bool high )
+{
+ if ( high == b_state )
+ return;
+
+ if (m_bUseLogicChain)
+ {
+ b_state = high;
+
+ for ( LogicIn * logic = this; logic; logic = logic->nextLogic() )
+ logic->setLastState(high);
+
+ if (m_bCanAddChanged)
+ {
+ m_pSimulator->addChangedLogic(this);
+ m_bCanAddChanged = false;
+ }
+
+ return;
+ }
+
+ m_old_g_out = m_g_out;
+ m_old_v_out = m_v_out;
+
+ // NOTE Make sure that the next two lines are the same as those in setElementSet and setLogic
+ m_g_out = high ? m_gHigh : m_gLow;
+ m_v_out = high ? m_vHigh : 0.0;
+
+ add_initial_dc();
+
+ m_old_g_out = 0.;
+ m_old_v_out = 0.;
+
+ b_state = high;
+
+ if ( p_eSet && p_eSet->circuit()->canAddChanged() )
+ {
+ m_pSimulator->addChangedCircuit( p_eSet->circuit() );
+ p_eSet->circuit()->setCanAddChanged(false);
+ }
+}
+
diff --git a/src/electronics/simulation/logic.h b/src/electronics/simulation/logic.h
new file mode 100644
index 0000000..be8374f
--- /dev/null
+++ b/src/electronics/simulation/logic.h
@@ -0,0 +1,211 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef LOGIC_H
+#define LOGIC_H
+
+#include "element.h"
+
+#include <qguardedptr.h>
+#include <qvaluelist.h>
+
+class Component;
+class Pin;
+class Simulator;
+
+typedef QValueList<QGuardedPtr<Pin> > PinList;
+
+typedef struct
+{
+ float risingTrigger; // Trigger on rising edge
+ float fallingTrigger; // Trigger on falling edge
+ float output; // Output voltage
+ float highImpedance; // Output impedance when high
+ float lowImpedance; // Output impedance when low
+} LogicConfig;
+
+class CallbackClass {};
+typedef void(CallbackClass::*CallbackPtr)( bool isHigh );
+
+/**
+Use this class for Logic Inputs - this will have infinite impedance.
+Use isHigh() will return whether the voltage level at the pin
+is high than the predetermined voltage threshold, and setHigh() will make the
+output high/low, also according to the predetermined logic type / voltages.
+
+@short Boolean Logic input
+*/
+class LogicIn : public Element
+{
+ public:
+ LogicIn( LogicConfig config );
+ virtual ~LogicIn();
+
+ virtual Type type() const { return Element_LogicIn; }
+ virtual void setElementSet( ElementSet *c );
+
+ /**
+ * Set logic values from the LogicConfig.
+ */
+ virtual void setLogic( LogicConfig config );
+ /**
+ * Check if the input state has changed, to see if we need to callback.
+ */
+ void check();
+ /**
+ * Returns whether the pin is 'high', as defined for the logic type
+ * Note: this is defined as the voltage on the pin, as opposed to what the
+ * state was set to (the two are not necessarily the same).
+ */
+ bool isHigh() const { return m_bLastState; }
+ /**
+ * When the logic state on this LogicIn changes, the function passed in this
+ * function will be called. At most one Callback can be added per LogicIn.
+ */
+ void setCallback( CallbackClass * object, CallbackPtr func );
+ /**
+ * Reads the LogicConfig values in from KTLConfig, and returns them in a
+ * nice object form.
+ */
+ static LogicConfig getConfig();
+ /**
+ * If this belongs to a logic chain, then this will be called from the chain.
+ */
+ void setLastState( bool state ) { m_bLastState = state; }
+ /**
+ * Returns a pointer to the next LogicIn in the chain.
+ */
+ LogicIn * nextLogic() const { return m_pNextLogic; }
+ /**
+ * Sets the next LogicIn in the chain.
+ */
+ void setNextLogic( LogicIn * next ) { m_pNextLogic = next; }
+ /**
+ * Calls the callback function, if there is one.
+ */
+ void callCallback()
+ {
+ if (m_pCallbackFunction)
+ (m_pCallbackObject->*m_pCallbackFunction)(m_bLastState);
+ }
+
+ protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+ CallbackPtr m_pCallbackFunction;
+ CallbackClass * m_pCallbackObject;
+ bool m_bLastState;
+ LogicIn * m_pNextLogic;
+ LogicConfig m_config;
+};
+
+
+/**
+@short Logic output/input
+*/
+class LogicOut : public LogicIn
+{
+ public:
+ LogicOut( LogicConfig config, bool _high );
+ virtual ~LogicOut();
+
+ virtual void setLogic( LogicConfig config );
+ virtual void setElementSet( ElementSet *c );
+ virtual void add_map();
+ virtual Type type() const { return Element_LogicOut; }
+
+ /**
+ * Call this function to override the logic-high output impedance as set by
+ * the user. Once set, the impedance will not be changed by the user
+ * updating the config; only by subsequent calls to this function.
+ */
+ void setOutputHighConductance( double g );
+ /**
+ * Call this function to override the logic-low output impedance as set by
+ * the user. Once set, the impedance will not be changed by the user
+ * updating the config; only by subsequent calls to this function.
+ */
+ void setOutputLowConductance( double g );
+ /**
+ * Call this function to override the logic out voltage as set by the
+ * user. Once set, the impedance will not be changed by the user
+ * updating the config; only by subsequent calls to this function.
+ */
+ void setOutputHighVoltage( double v );
+ /**
+ * Returns the voltage that this will output when high.
+ */
+ double outputHighVoltage() const { return m_vHigh; }
+ /**
+ * Sets the pin to be high/low
+ */
+ void setHigh( bool high );
+ /**
+ * @returns the state that this is outputting (regardless of voltage level on logic)
+ */
+ bool outputState() const { return b_state; }
+ /**
+ * Set whether or not this LogicOut is the head of a LogicChain (controls
+ * itself and a bunch of LogicIns).
+ */
+ void setUseLogicChain( bool use );
+ /**
+ * When a LogicOut configured as the start of a LogicChain changes start, it
+ * appends a pointer to itself to the list of change LogicOut, starting from
+ * the Simulator. This functions enables appending the next changed LogicOut
+ * to this one.
+ */
+ void setNextChanged( LogicOut * logicOut, unsigned char chain ) { m_pNextChanged[chain] = logicOut; }
+ /**
+ * To avoid a pointer to this LogicOut being added twice in one
+ * iteration due to the state changing twice, this LogicOut sets an
+ * added flag to true after adding it to the list of changed. The flag
+ * must be reset to false with this function (done by Simulator).
+ */
+ void setCanAddChanged( bool canAdd ) { m_bCanAddChanged = canAdd; }
+ /**
+ * Returns the next LogicOut that has changed, when configured as the start
+ * of a LogicChain.
+ * @see setNextChanged
+ */
+ LogicOut * nextChanged( unsigned char chain ) const { return m_pNextChanged[chain]; }
+ PinList pinList;
+ PinList::iterator pinListBegin;
+ PinList::iterator pinListEnd;
+
+ protected:
+ void configChanged();
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+ // Pre-initalized levels from config
+ double m_gHigh;
+ double m_gLow;
+ double m_vHigh;
+
+ // Whether to use the user-defined logic values
+ bool m_bOutputHighConductanceConst;
+ bool m_bOutputLowConductanceConst;
+ bool m_bOutputHighVoltageConst;
+
+ double m_g_out;
+ double m_v_out;
+ double m_old_g_out;
+ double m_old_v_out;
+ bool b_state;
+ bool m_bCanAddChanged;
+ LogicOut * m_pNextChanged[2];
+ Simulator * m_pSimulator;
+ bool m_bUseLogicChain;
+};
+
+#endif
+
diff --git a/src/electronics/simulation/matrix.cpp b/src/electronics/simulation/matrix.cpp
new file mode 100644
index 0000000..fb1248f
--- /dev/null
+++ b/src/electronics/simulation/matrix.cpp
@@ -0,0 +1,546 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "matrix.h"
+
+#include <kdebug.h>
+
+#include <assert.h>
+
+#include <cmath>
+#include <iostream>
+#include <vector>
+
+/// Minimum value before an entry is deemed "zero"
+const double epsilon = 1e-50;
+
+Matrix::Matrix( uint n, uint m )
+{
+ m_n = n;
+ m_m = m;
+ m_size = m_n+m_m;
+
+ m_mat = new matrix(m_size);
+ m_lu = new matrix(m_size);
+ m_y = new double[m_size];
+ m_inMap = new int[m_size];
+// m_outMap = new int[m_size];
+ m_map = new Map(m_size);
+ zero();
+}
+
+
+Matrix::~Matrix()
+{
+ delete m_map;
+ delete m_mat;
+ delete m_lu;
+ delete [] m_y;
+ delete [] m_inMap;
+// delete [] m_outMap;
+}
+
+
+void Matrix::zero()
+{
+ for ( uint i=0; i<m_size; i++ )
+ {
+ for ( uint j=0; j<m_size; j++ )
+ {
+ (*m_mat)[i][j] = 0.;
+ (*m_lu)[i][j] = 0.;
+ }
+ m_inMap[i] = i;
+// m_outMap[i] = i;
+ }
+
+ max_k = 0;
+}
+
+
+void Matrix::setUse( const uint i, const uint j, Map::e_type type, bool big )
+{
+ m_map->setUse( i, j, type, big );
+}
+
+
+void Matrix::createMap()
+{
+ int newMap[m_size];
+ m_map->createMap(newMap);
+ for ( uint i=0; i<m_size; i++ )
+ {
+ const int nu = newMap[i];
+ if ( nu != m_inMap[i] )
+ {
+ int old = -1;
+ for ( uint j=0; j<m_size && old == -1; j++ )
+ {
+ if ( m_inMap[j] == nu ) {
+ old = j;
+ }
+ }
+ assert( old != -1 );
+ swapRows( old, i );
+ }
+ }
+}
+
+
+void Matrix::swapRows( const uint a, const uint b )
+{
+ if ( a == b ) return;
+ m_mat->swapRows( a, b );
+
+ const int old = m_inMap[a];
+ m_inMap[a] = m_inMap[b];
+ m_inMap[b] = old;
+
+ max_k = 0;
+}
+
+
+/*void Matrix::genOutMap()
+{
+ for ( uint i=0; i<m_size; i++ )
+ {
+ m_outMap[ m_inMap[i] ] = i;
+ }
+}*/
+
+
+void Matrix::operator=( Matrix *const m )
+{
+ for ( uint _i=0; _i<m_size; _i++ )
+ {
+ uint i = m_inMap[_i];
+ for ( uint j=0; j<m_size; j++ )
+ {
+ (*m_mat)[i][j] = m->m(i,j);
+ }
+ }
+
+ max_k = 0;
+}
+
+void Matrix::operator+=( Matrix *const m )
+{
+ for ( uint _i=0; _i<m_size; _i++ )
+ {
+ uint i = m_inMap[_i];
+ for ( uint j=0; j<m_size; j++ )
+ {
+ (*m_mat)[i][j] += m->m(i,j);
+ }
+ }
+
+ max_k = 0;
+}
+
+void Matrix::performLU()
+{
+// max_k = 0;
+ uint n = m_size;
+ if ( n == 0 ) return;
+
+ // Copy the affected segment to LU
+ for ( uint i=max_k; i<n; i++ )
+ {
+ for ( uint j=max_k; j<n; j++ )
+ {
+ (*m_lu)[i][j] = (*m_mat)[i][j];
+ }
+ }
+
+ // LU decompose the matrix, and store result back in matrix
+ for ( uint k=0; k<n-1; k++ )
+ {
+ double * const lu_K_K = &(*m_lu)[k][k];
+ if ( std::abs(*lu_K_K) < 1e-10 )
+ {
+ if ( *lu_K_K < 0. ) *lu_K_K = -1e-10;
+ else *lu_K_K = 1e-10;
+ }
+ for ( uint i=std::max(k,max_k)+1; i<n; i++ )
+ {
+ (*m_lu)[i][k] /= *lu_K_K;
+ }
+ for ( uint i=std::max(k,max_k)+1; i<n; i++ )
+ {
+ const double lu_I_K = (*m_lu)[i][k];
+ if ( std::abs(lu_I_K) > 1e-12 )
+ {
+ for ( uint j=std::max(k,max_k)+1; j<n; j++ )
+ {
+ (*m_lu)[i][j] -= lu_I_K*(*m_lu)[k][j];
+ }
+ }
+ }
+ }
+
+ max_k = n;
+}
+
+void Matrix::fbSub( Vector* b )
+{
+ if ( m_size == 0 ) return;
+
+ for ( uint i=0; i<m_size; i++ )
+ {
+ m_y[m_inMap[i]] = (*b)[i];
+ }
+
+ // Forward substitution
+ for ( uint i=1; i<m_size; i++ )
+ {
+ double sum = 0;
+ for ( uint j=0; j<i; j++ )
+ {
+ sum += (*m_lu)[i][j]*m_y[j];
+ }
+ m_y[i] -= sum;
+ }
+
+ // Back substitution
+ m_y[m_size-1] /= (*m_lu)[m_size-1][m_size-1];
+ for ( int i=m_size-2; i>=0; i-- )
+ {
+ double sum = 0;
+ for ( uint j=i+1; j<m_size; j++ )
+ {
+ sum += (*m_lu)[i][j]*m_y[j];
+ }
+ m_y[i] -= sum;
+ m_y[i] /= (*m_lu)[i][i];
+ }
+
+ for ( uint i=0; i<m_size; i++ )
+ (*b)[i] = m_y[i];
+}
+
+
+void Matrix::multiply( Vector *rhs, Vector *result )
+{
+ if ( !rhs || !result ) return;
+ result->reset();
+ for ( uint _i=0; _i<m_size; _i++ )
+ {
+ uint i = m_inMap[_i];
+ for ( uint j=0; j<m_size; j++ )
+ {
+ (*result)[_i] += (*m_mat)[i][j] * (*rhs)[j];
+ }
+ }
+}
+
+
+void Matrix::displayMatrix()
+{
+ uint n = m_size;
+ for ( uint _i=0; _i<n; _i++ )
+ {
+ uint i = m_inMap[_i];
+ for ( uint j=0; j<n; j++ )
+ {
+ if ( j > 0 && (*m_mat)[i][j] >= 0 ) kdDebug() << "+";
+ kdDebug() << (*m_mat)[i][j] << "("<<j<<")";
+ }
+ kdDebug() << endl;
+ }
+}
+
+void Matrix::displayLU()
+{
+ uint n = m_size;
+ for ( uint _i=0; _i<n; _i++ )
+ {
+ uint i = m_inMap[_i];
+// uint i = _i;
+ for ( uint j=0; j<n; j++ )
+ {
+ if ( j > 0 && (*m_lu)[i][j] >= 0 ) std::cout << "+";
+ std::cout << (*m_lu)[i][j] << "("<<j<<")";
+ }
+ std::cout << std::endl;
+ }
+ std::cout << "m_inMap: ";
+ for ( uint i=0; i<n; i++ )
+ {
+ std::cout << i<<"->"<<m_inMap[i]<<" ";
+ }
+ std::cout << std::endl;
+ /*cout << "m_outMap: ";
+ for ( uint i=0; i<n; i++ )
+ {
+ cout << i<<"->"<<m_outMap[i]<<" ";
+ }
+ cout << endl;*/
+}
+
+
+Map::Map( const uint size )
+{
+ m_size = size;
+ m_map = new ETMap( m_size );
+ reset();
+}
+
+
+Map::~Map()
+{
+ delete m_map;
+}
+
+
+void Map::reset()
+{
+ for ( uint i=0; i<m_size; i++ )
+ {
+ for ( uint j=0; j<m_size; j++ )
+ {
+ (*m_map)[i][j] = 0;
+ }
+ }
+}
+
+
+void Map::setUse( const uint i, const uint j, Map::e_type type, bool big )
+{
+ if ( type == Map::et_none ) {
+ (*m_map)[i][j] = Map::et_none;
+ } else {
+ (*m_map)[i][j] = type | (big)?Map::et_big:0;
+ }
+}
+
+
+void Map::createMap( int *map )
+{
+ assert(map);
+
+ // In this function, the passes through that we make want to be done from
+ // top left to bottom right, to minimise fill-in
+
+ // available[i] is true if an external-row can be mapped to internal-row "i"
+ // map[i] gives the internal-row for external-row i
+ bool available[m_size];
+ for ( uint i=0; i<m_size; i++ )
+ {
+ available[i] = true;
+ map[i] = -1;
+ }
+
+ // This loop looks through columns and rows to find any swaps that are necessary
+ // (e.g. only one matrix-element in that row/column), and if no necessary swaps
+ // were found, then it will swap two rows according to criteria given below
+ bool badMap = false;
+ bool changed;
+ do
+ {
+ changed = false;
+
+ // Pass through columns
+ int E,N;
+ uint highest = 0;
+ for ( uint j=0; j<m_size; j++ )
+ {
+ if ( map[j] == -1 ) // If we haven't mapped this column yet
+ {
+ int count = 0; // Number of "spare" elements
+ int element; // Last element that is "spare", only applicable if count=1
+ for ( uint i=0; i<m_size; i++ )
+ {
+ if ( available[i] && (*m_map)[i][j] )
+ {
+ count++;
+ element = i;
+ }
+ }
+ if ( count == 0 ) {
+ badMap = true;
+ }
+ else if ( count == 1 )
+ {
+ const uint newType = (*m_map)[element][j];
+ if ( typeCmp( newType, highest) )
+ {
+ E=element;
+ N=j;
+ highest=newType;
+ }
+ }
+ }
+ }
+ // Pass through rows
+ for ( uint i=0; i<m_size; i++ )
+ {
+ if ( map[i] == -1 ) // If we haven't mapped this row yet
+ {
+ int count = 0; // Number of "spare" elements
+ int element; // Last element that is "spare", only applicable if count=1
+ for ( uint j=0; j<m_size; j++ )
+ {
+ if ( available[j] && (*m_map)[i][j] )
+ {
+ count++;
+ element = j;
+ }
+ }
+ if ( count == 0 ) {
+ badMap = true;
+ }
+ else if ( count == 1 )
+ {
+ const uint newType = (*m_map)[i][element];
+ if ( typeCmp( newType, highest) )
+ {
+ E=element;
+ N=i;
+ highest=newType;
+ }
+ }
+ }
+ }
+ if (highest)
+ {
+ available[E] = false;
+ map[N] = E;
+ changed = true;
+ }
+ if (!changed)
+ {
+ int next = -1; // next is the row to mapped to (interally)
+ uint j=0;
+
+ /// TODO We want to change this search to one that finds a swap, taking into acocunt the priorities given below
+ while ( next == -1 && j<m_size )
+ {
+ if ( available[j] ) next=j;
+ j++;
+ }
+ uint i=0;
+ while ( i<m_size && map[i] != -1 ) i++;
+ if ( next != -1 && i < m_size )
+ {
+ available[next] = false;
+ map[i] = next;
+ changed = true;
+ }
+ }
+ } while (changed);
+
+ if (badMap)
+ {
+// cerr << "Map::createMap: unable to create decent mapping; do not trust the matrix, Neo!"<<endl;
+ }
+
+ for ( int i = 0; i < int(m_size); ++i )
+ {
+ assert( map[i] >= 0 && map[i] < int(m_size) );
+ }
+
+ // Ignore this, for now:
+
+ // Now, we want to order the matrix, with the following priorities:
+ // (1) How often values change
+ // (2) How few values there are
+ // (3) How large the values are
+ // For each value in the column,
+}
+
+
+bool Map::typeCmp( const uint t1, const uint t2 )
+{
+ if (!t2) return true;
+ if (!t1) return false;
+
+ int t1_score = 1;
+ if ( t1 | Map::et_constant ) t1_score += 64;
+ else if ( t1 | Map::et_stable ) t1_score += 16;
+ else if ( t1 | Map::et_variable ) t1_score += 4;
+
+ int t2_score = 1;
+ if ( t2 | Map::et_constant ) t2_score += 64;
+ else if ( t2 | Map::et_stable ) t2_score += 16;
+ else if ( t2 | Map::et_variable ) t2_score += 4;
+
+ if ( t1 | Map::et_big ) t1_score *= 2;
+ if ( t2 | Map::et_big ) t2_score *= 2;
+
+ return ( t1_score >= t2_score );
+}
+
+
+Matrix22::Matrix22()
+{
+ reset();
+}
+
+bool Matrix22::solve()
+{
+ const double old_x1 = m_x1;
+ const double old_x2 = m_x2;
+
+ const bool e11 = std::abs((m_a11))<epsilon;
+ const bool e12 = std::abs((m_a12))<epsilon;
+ const bool e21 = std::abs((m_a21))<epsilon;
+ const bool e22 = std::abs((m_a22))<epsilon;
+
+ if (e11)
+ {
+ if ( e12||e21 )
+ return false;
+ m_x2 = m_b1/m_a12;
+ m_x1 = (m_b2-(m_a22*m_x2))/m_a21;
+ }
+ else if (e12)
+ {
+ if ( e11||e22 )
+ return false;
+ m_x1 = m_b1/m_a11;
+ m_x2 = (m_b2-(m_a21*m_x1))/m_a22;
+ }
+ else if (e21)
+ {
+ if ( e11||e22 )
+ return false;
+ m_x2 = m_b2/m_a22;
+ m_x1 = (m_b1-(m_a12*m_x2))/m_a11;
+ }
+ else if (e22)
+ {
+ if ( e12||e21 )
+ return false;
+ m_x1 = m_b2/m_a21;
+ m_x2 = (m_b1-(m_a11*m_x1))/m_a12;
+ }
+ else
+ {
+ m_x2 = (m_b2-(m_a21*m_b1/m_a11))/(m_a22-(m_a21*m_a12/m_a11));
+ m_x1 = (m_b1-(m_a12*m_x2))/m_a11;
+ }
+ if ( !std::isfinite(m_x1) || !std::isfinite(m_x2) )
+ {
+ m_x1 = old_x1;
+ m_x2 = old_x2;
+ return false;
+ }
+ return true;
+}
+
+void Matrix22::reset()
+{
+ m_a11=m_a12=m_a21=m_a22=0.;
+ m_b1=m_b2=0.;
+ m_x1=m_x2=0.;
+}
+
+
+
diff --git a/src/electronics/simulation/matrix.h b/src/electronics/simulation/matrix.h
new file mode 100644
index 0000000..4c3e518
--- /dev/null
+++ b/src/electronics/simulation/matrix.h
@@ -0,0 +1,248 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef MATRIX_H
+#define MATRIX_H
+
+#include "vec.h"
+
+#include <vector>
+
+class Map;
+
+typedef std::vector<std::vector<uint> > ETMap;
+
+/**
+@short Handles row-wise permutation of matrixes
+*/
+class Map
+{
+public:
+ enum e_type
+ {
+ et_none = 0x0, // Default value, none
+ et_constant = 0x1, // Never changes value in lifetime of matrix
+ et_stable = 0x2, // Value changes occasionally, e.g. user changing resistance value
+ et_variable = 0x3, // Rate of changing is unknown, probably average - e.g. logic
+ et_unstable = 0x4, // Value changes practically every call of performLU
+ et_big = 0x8 // Ignore this :-) It is used to OR with one of the others, but it should only ever be accessed by the class
+ };
+ Map( const uint size );
+ ~Map();
+
+ void setUse( const uint i, const uint j, Map::e_type type, bool big );
+ /**
+ * Generates an optimal permutation, returned in the given array
+ */
+ void createMap( int *map );
+ /**
+ * Resets the info to a blank pattern
+ */
+ void reset();
+
+protected:
+ /**
+ * Compares two "types", returning true if t1 >= t2, else false
+ */
+ bool typeCmp( const uint t1, const uint t2 );
+
+private:
+ ETMap *m_map;
+ uint m_size;
+};
+
+/**
+This class performs matrix storage, lu decomposition, forward and backward
+substitution, and a few other useful operations. Steps in using class:
+(1) Create an instance of this class with the correct size
+(2) Define the matrix pattern as neccessary:
+ (1) Call zero (unnecessary after initial ceration) to reset the pattern
+ & matrix
+ (2) Call setUse to set the use of each element in the matrix
+ (3) Call createMap to generate the row-wise permutation mapping for use
+ in partial pivoting
+(3) Add the values to the matrix
+(4) Call performLU, and get the results with fbSub
+(5) Repeat 2, 3, 4 or 5 as necessary.
+@todo We need to allow createMap to work while the matrix has already been initalised
+@short Matrix manipulation class tailored for circuit equations
+@author David Saxton
+*/
+class Matrix
+{
+public:
+ /**
+ * Creates a size x size square matrix m, with all values zero,
+ * and a right side vector x of size m+n
+ */
+ Matrix( uint n, uint m );
+ ~Matrix();
+ /**
+ * Sets all elements to zero
+ */
+ void zero();
+ /**
+ * Sets the type of (matrix-) element at the position i, j.
+ * @param type Describes how often the value changes
+ * @param big Set this true if the value is likely to be >= 1, else false
+ */
+ void setUse( const uint i, const uint j, Map::e_type type, bool big );
+ void setUse_b( const uint i, const uint j, Map::e_type type, bool big )
+ {
+ setUse( i, j+m_n, type, big );
+ }
+ void setUse_c( const uint i, const uint j, Map::e_type type, bool big )
+ {
+ setUse( i+m_n, j, type, big );
+ }
+ void setUse_d( const uint i, const uint j, Map::e_type type, bool big )
+ {
+ setUse( i+m_n, j+m_n, type, big );
+ }
+ /**
+ * Generates the row-wise permutation mapping from the values set by setUse
+ */
+ void createMap();
+ /**
+ * Returns true if the matrix is changed since last calling performLU()
+ * - i.e. if we do need to call performLU again.
+ */
+ inline bool isChanged() const { return max_k < m_size; }
+ /**
+ * Performs LU decomposition. Going along the rows,
+ * the value of the decomposed LU matrix depends only on
+ * the previous values.
+ */
+ void performLU();
+ /**
+ * Applies the right side vector (x) to the decomposed matrix,
+ * with the solution returned in x.
+ */
+ void fbSub( Vector* x );
+ /**
+ * Prints the matrix to stdout
+ */
+ void displayMatrix();
+ /**
+ * Prints the LU-decomposed matrix to stdout
+ */
+ void displayLU();
+ /**
+ * Sets the element matrix at row i, col j to value x
+ */
+ double& g( uint i, uint j )
+ {
+ i = m_inMap[i];
+ if ( i<max_k ) max_k=i;
+ if ( j<max_k ) max_k=j;
+
+ // I think I need the next line...
+ if ( max_k>0 ) max_k--;
+
+ return (*m_mat)[i][j];
+ }
+ double& b( uint i, uint j ) { return g( i, j+m_n ); }
+ double& c( uint i, uint j ) { return g( i+m_n, j ); }
+ double& d( uint i, uint j ) { return g( i+m_n, j+m_n ); }
+ const double g( uint i, uint j ) const { return (*m_mat)[m_inMap[i]][j]; }
+ const double b( uint i, uint j ) const { return g( i, j+m_n ); }
+ const double c( uint i, uint j ) const { return g( i+m_n, j ); }
+ const double d( uint i, uint j ) const { return g( i+m_n, j+m_n ); }
+ /**
+ * Returns the value of matrix at row i, col j.
+ */
+ const double m( uint i, uint j ) const
+ {
+ return (*m_mat)[m_inMap[i]][j];
+ }
+ /**
+ * Multiplies this matrix by the Vector rhs, and places the result
+ * in the vector pointed to by result. Will fail if wrong size.
+ */
+ void multiply( Vector *rhs, Vector *result );
+ /**
+ * Sets the values of this matrix to that of the given matrix
+ */
+ void operator=( Matrix *const m );
+ /**
+ * Adds the values of the given matrix to this matrix
+ */
+ void operator+=( Matrix *const m );
+
+private:
+ /**
+ * Swaps around the rows in the (a) the matrix; and (b) the mappings
+ */
+ void swapRows( const uint a, const uint b );
+
+// // Generates m_outMap from m_inMap
+// void genOutMap();
+
+ uint m_n;
+ uint m_m;
+ uint m_size;
+ uint max_k;
+
+ int *m_inMap; // Rowwise permutation mapping from external reference to internal storage
+// int *m_outMap; // Opposite of m_inMap
+
+ matrix *m_mat;
+ matrix *m_lu;
+ double *m_y; // Avoids recreating it lots of times
+ Map *m_map;
+};
+
+
+/**
+This class provides a very simple, lightweight, 2x2 matrix solver.
+It's fast and reliable. Set the values for the entries of A and b:
+
+A x = b
+
+call solve() (which returns true if successful - i.e. exactly one solution to the
+matrix), and get the values of x with the appropriate functions.
+
+@short 2x2 Matrix
+@author David Saxton
+*/
+class Matrix22
+{
+public:
+ Matrix22();
+
+ double &a11() { return m_a11; }
+ double &a12() { return m_a12; }
+ double &a21() { return m_a21; }
+ double &a22() { return m_a22; }
+
+ double &b1() { return m_b1; }
+ double &b2() { return m_b2; }
+
+ /**
+ * Solve the matrix. Returns true if successful (i.e. non-singular), else
+ * false. Get the solution with x1() and x2().
+ */
+ bool solve();
+ /**
+ * Resets all entries to zero
+ */
+ void reset();
+
+ double x1() const { return m_x1; }
+ double x2() const { return m_x2; }
+
+private:
+ double m_a11, m_a12, m_a21, m_a22;
+ double m_b1, m_b2;
+ double m_x1, m_x2;
+};
+
+
+#endif
diff --git a/src/electronics/simulation/nonlinear.cpp b/src/electronics/simulation/nonlinear.cpp
new file mode 100644
index 0000000..c975f16
--- /dev/null
+++ b/src/electronics/simulation/nonlinear.cpp
@@ -0,0 +1,97 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "matrix.h"
+#include "nonlinear.h"
+
+#include <cmath>
+using namespace std;
+
+const double KTL_MAX_DOUBLE = 1.7976931348623157e+308; ///< 7fefffff ffffffff
+const int KTL_MAX_EXPONENT = int( log( KTL_MAX_DOUBLE ) );
+
+
+NonLinear::NonLinear()
+ : Element()
+{
+}
+
+
+#ifndef MIN
+# define MIN(x,y) (((x) < (y)) ? (x) : (y))
+#endif
+
+
+// The function computes the exponential pn-junction current.
+double NonLinear::diodeCurrent( double v, double I_S, double Vte ) const
+{
+ return I_S * (exp( MIN( v / Vte, KTL_MAX_EXPONENT ) ) - 1);
+}
+
+
+
+double NonLinear::diodeConductance( double v, double I_S, double Vte ) const
+{
+ return I_S * exp( MIN( v / Vte, KTL_MAX_EXPONENT ) ) / Vte;
+}
+
+
+
+double NonLinear::diodeVoltage( double V, double V_prev, double V_T, double Vcrit ) const
+{
+ if ( V > Vcrit && fabs( V - V_prev ) > 2 * V_T )
+ {
+ if ( V_prev > 0 )
+ {
+ double arg = (V - V_prev) / V_T;
+ if (arg > 0)
+ V = V_prev + V_T * (2 + log( arg - 2 ));
+ else
+ V = V_prev - V_T * (2 + log( 2 - arg ));
+ }
+ else
+ V = (V_prev < 0) ? (V_T * log (V / V_T)) : Vcrit;
+ }
+ else
+ {
+ if ( V < 0 )
+ {
+ double arg = (V_prev > 0) ? (-1 - V_prev) : (2 * V_prev - 1);
+ if (V < arg)
+ V = arg;
+ }
+ }
+ return V;
+}
+
+
+double NonLinear::diodeCriticalVoltage( double I_S, double V_Te ) const
+{
+ return V_Te * log( V_Te / M_SQRT2 / I_S );
+}
+
+
+void NonLinear::diodeJunction( double V, double I_S, double V_Te, double * I, double * g ) const
+{
+ if (V < -3 * V_Te)
+ {
+ double a = 3 * V_Te / (V * M_E);
+ a = a * a * a;
+ *I = -I_S * (1 + a);
+ *g = +I_S * 3 * a / V;
+ }
+ else
+ {
+ double e = exp( MIN( V / V_Te, KTL_MAX_EXPONENT ) );
+ *I = I_S * (e - 1);
+ *g = I_S * e / V_Te;
+ }
+}
+
diff --git a/src/electronics/simulation/nonlinear.h b/src/electronics/simulation/nonlinear.h
new file mode 100644
index 0000000..422f3a4
--- /dev/null
+++ b/src/electronics/simulation/nonlinear.h
@@ -0,0 +1,48 @@
+/***************************************************************************
+ * Copyright (C) 2003-2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef NONLINEAR_H
+#define NONLINEAR_H
+
+#include "element.h"
+
+/**
+@short Represents a non-linear circuit element (such as a diode)
+@author David Saxton
+*/
+class NonLinear : public Element
+{
+ public:
+ NonLinear();
+
+ virtual bool isNonLinear() { return true; }
+ /**
+ * Newton-Raphson iteration: Update equation system.
+ */
+ virtual void update_dc() = 0;
+
+ protected:
+ double diodeCurrent( double v, double I_S, double Vte ) const;
+ /**
+ * Conductance of the diode - the derivative of Schockley's
+ * approximation.
+ */
+ double diodeConductance( double v, double I_S, double Vte ) const;
+ /**
+ * Limits the diode voltage to prevent divergence in the nonlinear
+ * iterations.
+ */
+ double diodeVoltage( double v, double V_prev, double Vt, double V_crit ) const;
+ void diodeJunction( double v, double I_S, double Vte, double * I, double * g ) const;
+
+ double diodeCriticalVoltage( double I_S, double Vte ) const;
+};
+
+#endif
diff --git a/src/electronics/simulation/opamp.cpp b/src/electronics/simulation/opamp.cpp
new file mode 100644
index 0000000..45fbf02
--- /dev/null
+++ b/src/electronics/simulation/opamp.cpp
@@ -0,0 +1,77 @@
+/***************************************************************************
+ * Copyright (C) 2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementset.h"
+#include "matrix.h"
+#include "opamp.h"
+
+OpAmp::OpAmp()
+ : Element::Element()
+{
+ m_numCBranches = 1;
+ m_numCNodes = 3;
+}
+
+
+OpAmp::~OpAmp()
+{
+}
+
+
+void OpAmp::add_map()
+{
+ if (!b_status)
+ return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ // Non-inverting input
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_constant, true );
+ }
+
+ if ( !p_cnode[2]->isGround )
+ {
+ // Inverting input
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[2]->n(), Map::et_constant, true );
+ }
+
+ if ( !p_cnode[1]->isGround )
+ {
+ // Output
+ p_A->setUse_b( p_cnode[1]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ }
+}
+
+
+void OpAmp::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ // Non-inverting input
+ A_c( 0, 0 ) = 1;
+
+ // Inverting input
+ A_c( 0, 2 ) = -1;
+
+ // Output
+ A_b( 1, 0 ) = 1;
+}
+
+
+void OpAmp::updateCurrents()
+{
+ if (!b_status) return;
+ m_cnodeI[0] = m_cnodeI[2] = 0.0;
+ m_cnodeI[1] = p_cbranch[0]->i;
+}
+
+
+
diff --git a/src/electronics/simulation/opamp.h b/src/electronics/simulation/opamp.h
new file mode 100644
index 0000000..22fe843
--- /dev/null
+++ b/src/electronics/simulation/opamp.h
@@ -0,0 +1,36 @@
+/***************************************************************************
+ * Copyright (C) 2005 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef OPAMP_H
+#define OPAMP_H
+
+#include <element.h>
+
+/**
+node 0: non-inverting input
+node 1: output
+node 2: inverting input
+@author David Saxton
+*/
+class OpAmp : public Element
+{
+ public:
+ OpAmp();
+ virtual ~OpAmp();
+
+ virtual Type type() const { return Element_OpAmp; }
+ virtual void add_map();
+
+ protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+};
+
+#endif
diff --git a/src/electronics/simulation/reactive.cpp b/src/electronics/simulation/reactive.cpp
new file mode 100644
index 0000000..83fcfd4
--- /dev/null
+++ b/src/electronics/simulation/reactive.cpp
@@ -0,0 +1,36 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+
+#include "reactive.h"
+
+Reactive::Reactive( const double delta )
+ : Element()
+{
+ m_delta = delta;
+}
+
+
+Reactive::~Reactive()
+{
+}
+
+
+void Reactive::setDelta( double delta )
+{
+ m_delta = delta;
+ updateStatus();
+}
+
+
+bool Reactive::updateStatus()
+{
+ return Element::updateStatus();
+}
diff --git a/src/electronics/simulation/reactive.h b/src/electronics/simulation/reactive.h
new file mode 100644
index 0000000..1142b34
--- /dev/null
+++ b/src/electronics/simulation/reactive.h
@@ -0,0 +1,42 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef REACTIVE_H
+#define REACTIVE_H
+
+#include "element.h"
+
+/**
+@short Represents a reactive element (such as a capacitor)
+@author David Saxton
+*/
+class Reactive : public Element
+{
+public:
+ Reactive( const double delta );
+ virtual ~Reactive();
+
+ virtual bool isReactive() { return true; }
+ /**
+ * Call this function to set the time period (in seconds)
+ */
+ void setDelta( double delta );
+ /**
+ * Called on every time step for the element to update itself
+ */
+ virtual void time_step() = 0;
+
+protected:
+ virtual bool updateStatus();
+
+ double m_delta; // Delta time interval
+};
+
+#endif
diff --git a/src/electronics/simulation/resistance.cpp b/src/electronics/simulation/resistance.cpp
new file mode 100644
index 0000000..9c2d25d
--- /dev/null
+++ b/src/electronics/simulation/resistance.cpp
@@ -0,0 +1,95 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementset.h"
+#include "matrix.h"
+#include "resistance.h"
+
+// #include <kdebug.h>
+
+Resistance::Resistance( const double resistance )
+ : Element::Element()
+{
+ m_g = resistance < 1e-9 ? 1e9 : 1./resistance;
+ m_numCNodes = 2;
+// kdDebug() << k_funcinfo << endl;
+}
+
+Resistance::~Resistance()
+{
+// kdDebug() << k_funcinfo << endl;
+}
+
+void Resistance::setConductance( const double g )
+{
+ if ( g == m_g )
+ return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ // Remove old resistance
+ m_g = -m_g;
+ add_initial_dc();
+
+ m_g = g;
+ add_initial_dc();
+}
+
+
+void Resistance::setResistance( const double r )
+{
+ setConductance( r < 1e-9 ? 1e9 : 1./r );
+}
+
+
+void Resistance::add_map()
+{
+ if (!b_status)
+ return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[0]->n(), Map::et_stable, false );
+ }
+ if ( !p_cnode[1]->isGround ) {
+ p_A->setUse( p_cnode[1]->n(), p_cnode[1]->n(), Map::et_stable, false );
+ }
+
+ if ( !p_cnode[0]->isGround && !p_cnode[1]->isGround )
+ {
+ p_A->setUse( p_cnode[0]->n(), p_cnode[1]->n(), Map::et_stable, false );
+ p_A->setUse( p_cnode[1]->n(), p_cnode[0]->n(), Map::et_stable, false );
+ }
+}
+
+
+void Resistance::add_initial_dc()
+{
+ if (!b_status) return;
+
+ A_g( 0, 0 ) += m_g;
+ A_g( 1, 1 ) += m_g;
+ A_g( 0, 1 ) -= m_g;
+ A_g( 1, 0 ) -= m_g;
+}
+
+
+void Resistance::updateCurrents()
+{
+ if (!b_status) return;
+ const double v=p_cnode[0]->v-p_cnode[1]->v;
+ m_cnodeI[1] = v*m_g;
+ m_cnodeI[0] = -m_cnodeI[1];
+}
+
+
+
+
diff --git a/src/electronics/simulation/resistance.h b/src/electronics/simulation/resistance.h
new file mode 100644
index 0000000..7e5a62e
--- /dev/null
+++ b/src/electronics/simulation/resistance.h
@@ -0,0 +1,43 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef RESISTANCE_H
+#define RESISTANCE_H
+
+#include "element.h"
+
+/**
+@short Resistance
+@author David saxton
+*/
+class Resistance : public Element
+{
+public:
+ Resistance( const double resistance );
+ virtual ~Resistance();
+
+ virtual Type type() const { return Element_Resistance; }
+
+ void setConductance( const double g );
+ void setResistance( const double r );
+
+ double resistance() { return 1/m_g; }
+ double conductance() { return m_g; }
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_g; // Conductance
+};
+
+#endif
diff --git a/src/electronics/simulation/vccs.cpp b/src/electronics/simulation/vccs.cpp
new file mode 100644
index 0000000..7ff1bc1
--- /dev/null
+++ b/src/electronics/simulation/vccs.cpp
@@ -0,0 +1,93 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementset.h"
+#include "matrix.h"
+#include "vccs.h"
+
+VCCS::VCCS( const double gain )
+ : Element::Element()
+{
+ m_g = gain;
+ m_numCBranches = 1;
+ m_numCNodes = 4;
+}
+
+
+VCCS::~VCCS()
+{
+}
+
+
+void VCCS::setGain( const double g )
+{
+ if ( g == m_g ) return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ // Remove old values
+ m_g = -m_g;
+ add_initial_dc();
+
+ // Add new values
+ m_g = g;
+ add_initial_dc();
+}
+
+
+void VCCS::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_stable, false );
+ }
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[1]->n(), Map::et_stable, false );
+ }
+ if ( !p_cnode[2]->isGround )
+ {
+ p_A->setUse_b( p_cnode[2]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[2]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[3]->isGround )
+ {
+ p_A->setUse_b( p_cnode[3]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[3]->n(), Map::et_constant, true );
+ }
+}
+
+
+void VCCS::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ A_g( 2, 0 ) += m_g;
+ A_g( 3, 0 ) -= m_g;
+ A_g( 2, 1 ) -= m_g;
+ A_g( 3, 1 ) += m_g;
+}
+
+
+void VCCS::updateCurrents()
+{
+ if (!b_status)
+ return;
+
+ m_cnodeI[0] = m_cnodeI[1] = 0.;
+ m_cnodeI[3] = (p_cnode[0]->v-p_cnode[1]->v)*m_g;
+ m_cnodeI[2] = -m_cnodeI[3];
+}
+
+
diff --git a/src/electronics/simulation/vccs.h b/src/electronics/simulation/vccs.h
new file mode 100644
index 0000000..26f1101
--- /dev/null
+++ b/src/electronics/simulation/vccs.h
@@ -0,0 +1,40 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef VCCS_H
+#define VCCS_H
+
+#include "element.h"
+
+/**
+CNodes n0 and n1 are used for the voltage control.
+CNodes n2 and n3 are used for the current output.
+@short Voltage Controlled Current Source
+@author David Saxton
+*/
+class VCCS : public Element
+{
+public:
+ VCCS( const double gain );
+ virtual ~VCCS();
+
+ virtual Type type() const { return Element_VCCS; }
+ void setGain( const double g );
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_g; // Conductance
+};
+
+#endif
diff --git a/src/electronics/simulation/vcvs.cpp b/src/electronics/simulation/vcvs.cpp
new file mode 100644
index 0000000..68604df
--- /dev/null
+++ b/src/electronics/simulation/vcvs.cpp
@@ -0,0 +1,93 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "elementset.h"
+#include "matrix.h"
+#include "vcvs.h"
+
+VCVS::VCVS( const double gain )
+ : Element::Element()
+{
+ m_g = gain;
+ m_numCBranches = 1;
+ m_numCNodes = 4;
+}
+
+
+VCVS::~VCVS()
+{
+}
+
+
+void VCVS::setGain( const double g )
+{
+ if ( g == m_g )
+ return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ m_g = -m_g;
+ add_initial_dc();
+
+ m_g = g;
+ add_initial_dc();
+}
+
+
+void VCVS::add_map()
+{
+ if (!b_status)
+ return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_stable, false );
+ }
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[1]->n(), Map::et_stable, false );
+ }
+ if ( !p_cnode[2]->isGround )
+ {
+ p_A->setUse_b( p_cnode[2]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[2]->n(), Map::et_constant, true );
+ }
+ if ( !p_cnode[3]->isGround )
+ {
+ p_A->setUse_b( p_cnode[3]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[3]->n(), Map::et_constant, true );
+ }
+}
+
+
+void VCVS::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ A_c( 0, 0 ) -= m_g;
+ A_c( 0, 1 ) += m_g;
+ A_b( 2, 0 ) = 1;
+ A_c( 0, 2 ) = 1;
+ A_b( 3, 0 ) = -1;
+ A_c( 0, 3 ) = -1;
+}
+
+
+void VCVS::updateCurrents()
+{
+ if (!b_status) return;
+ m_cnodeI[0] = m_cnodeI[1] = 0.;
+ m_cnodeI[3] = p_cbranch[0]->i;
+ m_cnodeI[2] = -m_cnodeI[3];
+}
+
+
diff --git a/src/electronics/simulation/vcvs.h b/src/electronics/simulation/vcvs.h
new file mode 100644
index 0000000..862dea9
--- /dev/null
+++ b/src/electronics/simulation/vcvs.h
@@ -0,0 +1,40 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef VCVS_H
+#define VCVS_H
+
+#include "element.h"
+
+/**
+Voltage source between nodes c2 and c3
+Controlling voltage between nodes c0 and c1
+@short Voltage Controlled Voltage Source
+@author David Saxton
+*/
+class VCVS : public Element
+{
+public:
+ VCVS( const double gain );
+ virtual ~VCVS();
+
+ virtual Type type() const { return Element_VCVS; }
+ void setGain( const double g );
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_g; // Conductance
+};
+
+#endif
diff --git a/src/electronics/simulation/vec.cpp b/src/electronics/simulation/vec.cpp
new file mode 100644
index 0000000..d45f505
--- /dev/null
+++ b/src/electronics/simulation/vec.cpp
@@ -0,0 +1,166 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "vec.h"
+
+#include <assert.h>
+#include <cmath>
+#include <string.h>
+using namespace std;
+
+Vector::Vector( const int size )
+{
+ m_size = size;
+ m_vec = new double[m_size];
+ reset();
+ b_changed = true;
+}
+
+
+Vector::~Vector()
+{
+ // Hmm...this looks like it's the correct format, although valgrind complains
+ // about improper memory use. Interesting. "delete m_vec" definitely leaks
+ // memory, so this seems like the lesser of two evils. I miss C memory allocation
+ // somtimes, with a nice simple free :p
+ delete [] m_vec;
+}
+
+
+void Vector::reset()
+{
+ for ( int i=0; i<m_size; i++ )
+ {
+ m_vec[i] = 0.;
+ }
+ b_changed = true;
+}
+
+
+void Vector::limitTo( double scaleMax, Vector * limitVector )
+{
+ assert( limitVector );
+ assert( limitVector->size() == size() );
+
+ for ( int i = 0; i < m_size; ++i )
+ {
+ double limitAbs = std::abs( limitVector->m_vec[i] );
+ if ( limitAbs < 1e-6 )
+ limitAbs = 1e-6;
+
+ double thisAbs = std::abs( m_vec[i] );
+ if ( thisAbs < 1e-6 )
+ thisAbs = 1e-6;
+
+ if ( (thisAbs / limitAbs) > scaleMax )
+ m_vec[i] /= (thisAbs / limitAbs);
+
+ else if ( (limitAbs / thisAbs) > scaleMax )
+ m_vec[i] /= (limitAbs / thisAbs);
+ }
+ b_changed = true;
+}
+
+
+void Vector::operator += ( Vector *rhs )
+{
+ if (!rhs) return;
+ for ( int i=0; i<m_size; i++ )
+ {
+ m_vec[i] += (*rhs)[i];
+ }
+ b_changed = true;
+}
+
+
+void Vector::operator -= ( Vector *rhs )
+{
+ if (!rhs) return;
+ for ( int i=0; i<m_size; i++ )
+ {
+ m_vec[i] -= (*rhs)[i];
+ }
+ b_changed = true;
+}
+
+
+void Vector::operator *=( double s )
+{
+ for ( int i=0; i<m_size; i++ )
+ {
+ m_vec[i] *= s;
+ }
+ b_changed = true;
+}
+
+
+void Vector::operator = ( Vector& v )
+{
+ assert( size() == v.size() );
+ memcpy( m_vec, v.m_vec, m_size * sizeof( double ) );
+ b_changed = true;
+}
+
+
+void Vector::negative( Vector *rhs )
+{
+ if (!rhs) return;
+ for ( int i=0; i<m_size; i++ )
+ {
+ m_vec[i] = -(*rhs)[i];
+ }
+ b_changed = true;
+}
+
+
+double Vector::abs() const
+{
+ double s=0;
+ for ( int i=0; i<m_size; i++ )
+ {
+ s += m_vec[i]*m_vec[i];
+ }
+ return sqrt(s);
+}
+
+
+
+// matrix stuff...
+
+matrix::matrix( const uint size )
+{
+ m_size = size;
+ m_mat = new Vector*[m_size];
+ for ( uint i=0; i<m_size; ++i )
+ {
+ m_mat[i] = new Vector(m_size);
+ }
+}
+
+matrix::~matrix()
+{
+ for ( uint i=0; i<m_size; ++i )
+ {
+ delete m_mat[i];
+ }
+ delete [] m_mat;
+}
+
+
+void matrix::swapRows( const uint a, const uint b )
+{
+ if ( a == b ) return;
+ Vector *v = m_mat[a];
+ m_mat[a] = m_mat[b];
+ m_mat[b] = v;
+}
+
+
+
diff --git a/src/electronics/simulation/vec.h b/src/electronics/simulation/vec.h
new file mode 100644
index 0000000..7192bb3
--- /dev/null
+++ b/src/electronics/simulation/vec.h
@@ -0,0 +1,109 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef VEC_H
+#define VEC_H
+
+typedef unsigned uint;
+
+/**
+@short Vector of doubles, faster than STL vector
+@author David Saxton
+*/
+class Vector
+{
+public:
+ Vector( const int size );
+ ~Vector();
+
+ double & operator[]( const int i ) { b_changed=true; return m_vec[i]; }
+ const double operator[]( const int i ) const { return m_vec[i]; }
+ int size() const { return m_size; }
+ /**
+ * Resets all values to 0
+ */
+ void reset();
+ /**
+ * Limits the absolute value of each component of this vector to scaleMax
+ * times bigger or smaller than the components of limitVector.
+ */
+ void limitTo( double scaleMax, Vector * limitVector );
+ /**
+ * Adds the Vector rhs to this
+ */
+ void operator+=( Vector *rhs );
+ /**
+ * Subtracts the Vector rhs from this
+ */
+ void operator-=( Vector *rhs );
+ /**
+ * Multiplies this Vector by the given scaler constant
+ */
+ void operator*=( double s );
+ /**
+ * Sets this vector equal to the given vector
+ */
+ void operator=( Vector& v );
+ /**
+ * Copies the negative values of the given vector to this vector.
+ * (i.e. sets this = -rhs )
+ */
+ void negative( Vector *rhs );
+ /**
+ * Returns the absolute value of this vector, defined as the squareroot
+ * of the sum of the square of each element of the vector
+ */
+ double abs() const;
+ /**
+ * Returns true if the vector has changed since setUnchanged was last called
+ * Note that this will return true if the vector has just been read, due to
+ * limitations with the [] operator.
+ */
+ inline bool isChanged() const { return b_changed; }
+ /**
+ * Sets the changed status to false.
+ */
+ inline void setUnchanged() { b_changed=false; }
+
+private:
+ Vector( const Vector & );
+ Vector & operator= ( const Vector & );
+
+ bool b_changed;
+ double *m_vec;
+ int m_size;
+};
+
+/**
+@short Container for Vector of Vectors
+@author David Saxton
+*/
+class matrix
+{
+public:
+ matrix( const uint size );
+ ~matrix();
+
+ Vector & operator[]( const uint i ) { return *(m_mat[i]); }
+ const Vector & operator[]( const uint i ) const { return *(m_mat[i]); }
+ /**
+ * Swaps the pointers to the given rows
+ */
+ void swapRows( const uint a, const uint b );
+
+private:
+ matrix( const matrix & );
+ matrix & operator= ( const matrix & );
+
+ Vector **m_mat;
+ uint m_size;
+};
+
+#endif
diff --git a/src/electronics/simulation/voltagepoint.cpp b/src/electronics/simulation/voltagepoint.cpp
new file mode 100644
index 0000000..44d3b01
--- /dev/null
+++ b/src/electronics/simulation/voltagepoint.cpp
@@ -0,0 +1,73 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "matrix.h"
+#include "voltagepoint.h"
+#include "elementset.h"
+
+VoltagePoint::VoltagePoint( const double voltage )
+ : Element::Element()
+{
+ m_voltage = -voltage;
+ m_numCBranches = 1;
+ m_numCNodes = 1;
+}
+
+
+VoltagePoint::~VoltagePoint()
+{
+}
+
+
+void VoltagePoint::setVoltage( const double v )
+{
+ if ( -v == m_voltage ) return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ m_voltage = -v;
+ add_initial_dc();
+}
+
+
+void VoltagePoint::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_b( p_cnode[0]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_constant, true );
+ }
+}
+
+
+void VoltagePoint::add_initial_dc()
+{
+ if (!b_status) return;
+
+ A_b( 0, 0 ) = -1;
+ A_c( 0, 0 ) = -1;
+
+ b_v( 0 ) = m_voltage;
+}
+
+
+void VoltagePoint::updateCurrents()
+{
+ if (!b_status) return;
+ m_cnodeI[0] = p_cbranch[0]->i;
+}
+
+
+
+
+
diff --git a/src/electronics/simulation/voltagepoint.h b/src/electronics/simulation/voltagepoint.h
new file mode 100644
index 0000000..924b0af
--- /dev/null
+++ b/src/electronics/simulation/voltagepoint.h
@@ -0,0 +1,39 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef VOLTAGEPOINT_H
+#define VOLTAGEPOINT_H
+
+#include "element.h"
+
+/**
+@short VoltagePoint
+@author David saxton
+*/
+class VoltagePoint : public Element
+{
+public:
+ VoltagePoint( const double voltage );
+ virtual ~VoltagePoint();
+
+ virtual Type type() const { return Element_VoltagePoint; }
+ void setVoltage( const double voltage );
+ double voltage() { return m_voltage; }
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_voltage; // Conductance
+};
+
+#endif
diff --git a/src/electronics/simulation/voltagesignal.cpp b/src/electronics/simulation/voltagesignal.cpp
new file mode 100644
index 0000000..d7a5839
--- /dev/null
+++ b/src/electronics/simulation/voltagesignal.cpp
@@ -0,0 +1,79 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "matrix.h"
+#include "voltagesignal.h"
+#include "elementset.h"
+
+VoltageSignal::VoltageSignal( const double delta, const double voltage )
+ : Reactive::Reactive(delta)
+{
+ m_voltage = voltage;
+ m_numCNodes = 2;
+ m_numCBranches = 1;
+}
+
+
+VoltageSignal::~VoltageSignal()
+{
+}
+
+
+void VoltageSignal::setVoltage( const double v )
+{
+ m_voltage = v;
+}
+
+
+void VoltageSignal::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_b( p_cnode[0]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_constant, true );
+ }
+
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse_b( p_cnode[1]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[1]->n(), Map::et_constant, true );
+ }
+}
+
+
+void VoltageSignal::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ A_b( 0, 0 ) = -1;
+ A_c( 0, 0 ) = -1;
+ A_b( 1, 0 ) = 1;
+ A_c( 0, 1 ) = 1;
+}
+
+
+void VoltageSignal::time_step()
+{
+ if (!b_status) return;
+ b_v( 0 ) = m_voltage*advance();
+}
+
+
+void VoltageSignal::updateCurrents()
+{
+ if (!b_status) return;
+ m_cnodeI[1] = p_cbranch[0]->i;
+ m_cnodeI[0] = -m_cnodeI[1];
+}
+
+
diff --git a/src/electronics/simulation/voltagesignal.h b/src/electronics/simulation/voltagesignal.h
new file mode 100644
index 0000000..5e35e72
--- /dev/null
+++ b/src/electronics/simulation/voltagesignal.h
@@ -0,0 +1,41 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef VOLTAGESIGNAL_H
+#define VOLTAGESIGNAL_H
+
+#include "reactive.h"
+#include "elementsignal.h"
+
+/**
+@short VoltageSignal
+@author David saxton
+*/
+class VoltageSignal : public Reactive, public ElementSignal
+{
+public:
+ VoltageSignal( const double delta, const double voltage );
+ virtual ~VoltageSignal();
+
+ virtual Element::Type type() const { return Element_VoltageSignal; }
+ void setVoltage( const double voltage );
+ double voltage() { return m_voltage; }
+ virtual void time_step();
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_voltage; // Voltage
+};
+
+#endif
diff --git a/src/electronics/simulation/voltagesource.cpp b/src/electronics/simulation/voltagesource.cpp
new file mode 100644
index 0000000..b4fca64
--- /dev/null
+++ b/src/electronics/simulation/voltagesource.cpp
@@ -0,0 +1,77 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#include "matrix.h"
+#include "voltagesource.h"
+#include "elementset.h"
+
+VoltageSource::VoltageSource( const double voltage )
+ : Element::Element()
+{
+ m_v = voltage;
+ m_numCBranches = 1;
+ m_numCNodes = 2;
+}
+
+VoltageSource::~VoltageSource()
+{
+}
+
+void VoltageSource::setVoltage( const double v )
+{
+ if ( m_v == v ) return;
+
+ if (p_eSet)
+ p_eSet->setCacheInvalidated();
+
+ m_v = v;
+ add_initial_dc();
+}
+
+
+void VoltageSource::add_map()
+{
+ if (!b_status) return;
+
+ if ( !p_cnode[0]->isGround )
+ {
+ p_A->setUse_b( p_cnode[0]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[0]->n(), Map::et_constant, true );
+ }
+
+ if ( !p_cnode[1]->isGround )
+ {
+ p_A->setUse_b( p_cnode[1]->n(), p_cbranch[0]->n(), Map::et_constant, true );
+ p_A->setUse_c( p_cbranch[0]->n(), p_cnode[1]->n(), Map::et_constant, true );
+ }
+}
+
+
+void VoltageSource::add_initial_dc()
+{
+ if (!b_status)
+ return;
+
+ A_b( 0, 0 ) = -1;
+ A_c( 0, 0 ) = -1;
+ A_b( 1, 0 ) = 1;
+ A_c( 0, 1 ) = 1;
+
+ b_v( 0 ) = m_v;
+}
+
+void VoltageSource::updateCurrents()
+{
+ if (!b_status) return;
+ m_cnodeI[0] = p_cbranch[0]->i;
+ m_cnodeI[1] = -m_cnodeI[0];
+}
+
+
diff --git a/src/electronics/simulation/voltagesource.h b/src/electronics/simulation/voltagesource.h
new file mode 100644
index 0000000..9b27b0d
--- /dev/null
+++ b/src/electronics/simulation/voltagesource.h
@@ -0,0 +1,38 @@
+/***************************************************************************
+ * Copyright (C) 2003-2004 by David Saxton *
+ * david@bluehaze.org *
+ * *
+ * This program is free software; you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation; either version 2 of the License, or *
+ * (at your option) any later version. *
+ ***************************************************************************/
+
+#ifndef VOLTAGESOURCE_H
+#define VOLTAGESOURCE_H
+
+#include "element.h"
+
+/**
+CNode n0 is the negative terminal, CNode n1 is the positive terminal
+@short Voltage Source
+*/
+class VoltageSource : public Element
+{
+public:
+ VoltageSource( const double voltage );
+ virtual ~VoltageSource();
+
+ virtual Type type() const { return Element_VoltageSource; }
+ void setVoltage( const double v );
+ virtual void add_map();
+
+protected:
+ virtual void updateCurrents();
+ virtual void add_initial_dc();
+
+private:
+ double m_v; // Voltage
+};
+
+#endif