diff options
Diffstat (limited to 'src/electronics/simulation')
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 |