Ocs package: Difference between revisions

Jump to navigation Jump to search

Ocs package (edit)

Revision as of 19:02, 2 September 2018

12,108 bytes added ,  2 September 2018
m
Rename "Octave-Forge" to "Octave Forge" (https://lists.gnu.org/archive/html/octave-maintainers/2018-08/msg00138.html).
m (Rename "Octave-Forge" to "Octave Forge" (https://lists.gnu.org/archive/html/octave-maintainers/2018-08/msg00138.html).)
(33 intermediate revisions by 3 users not shown)
Line 2: Line 2:
__TOC__
__TOC__
== History and Motivation ==
== History and Motivation ==
OCS was developed during the CoMSON (Coupled Multiscale Simulation and Optimization) project which involved several universities but also several industrial partners.
Each of the industrial partners at the time was using its own circuit simulation software and each software had different file formats for circuit netlists.
Given the purposes of the project and the composition of the consortium the main design objectives for OCS where
* provide a format for "element evaluators" independent of time-stepping algorithms
* provide a "hierarchical" data structure where elements could be composed themselves of lumped-element networks
* allow coupling of lumped-element networks (0D) and 1D/2D/3D device models
* use an intermediate/interchange file format so that none of the formats in use by the industrial partners would be favoured over the others
* be written in an interpreted language for quick prototyping and easy maintainance
* be Free Software
== Problem Formulation ==
== Problem Formulation ==
The circuit description in OCS is based on (a variant of) modified nodal analysis (MNA) model for lumped-element networks.
It is easy to verify that the common charge/flux-based MNA model is a special case of the model presented below.
We consider a circuit with M elements and N nodes, the core of the MNA model is a set of N equations of the form
<math>
\sum_{{m}=1}^{M}  F_{mn} =  0
\qquad
n = 1, \, \ldots \, ,N
</math>
where <math>F_{mn}</math> denotes the current from the node n due to element m.
The equations above are the Kirchhoff current law (KCL) for each of the electrical nodes of the network.
The currents can be expressed in terms of the node voltages <math>e</math> and the internal variables <math>r_m \; (m = 1\ldots M)</math>
<math>
F_{mn} =
A_{mn} \dot{r}_{m} + J_{mn} \left({e}, {r}_{m} \right)
\qquad
\begin{array}{l}
n = 1, \, \ldots \, ,N \\
m = 1, \, \ldots \, ,M
\end{array}
</math>
Notice that the variables <math>{{r}}_{m}</math> only appear in the equations defining the fluxes relative to the m-th element, for this reason they are sometimes referred to as internal variables of the m-th element.
The full MNA model is finally obtained by substituting the current definitions in the KCL and complementing it with a suitable number <math>I_{m}</math> of constitutive relations for the internal variables of each element
<math>
\sum_{{m}=1}^{M} \left[
\ A_{mn} \dot{{r}}_{m} +
J_{mn} \left( {e}, {r}_{m} \right)
\right] = 0
\qquad
\begin{array}{l}
{n} = 1, \, \ldots \, ,N
\end{array}
</math>
<math>
B_{mi} \dot{{r}}_{m} +
Q_{mi}\ \left( {e}, {r}_{m}  \right) = 0 \qquad
\begin{array}{l}
{i} = 1, \, \ldots \, ,{I}_m \\
{m} = 1, \, \ldots \, ,M
\end{array}
</math>
Notice that the assumption that only time derivatives of internal variables appear above and that terms involving such derivatives are linear does not impose restrictions on the applicability of the model.
== Data Structure ==
== Data Structure ==
A circuit is represented in OCS by a struct variable with the fields listed below
{{Code|OCS structure format |<syntaxhighlight lang="octave" style="font-size:13px">
cir_struct =
{
  LCR:  struct      % the fields of LCR are shown below
  NLC:  struct      % NLC has the same fields as LCR
  namesn: matrix    % numbers of vars that are assigned a name in and.nms
  namess: cell      % the names corresponding to the vars above
  totextvar: scalar % the total number of external variables
  totintvar: scalar % the total number of internal variables
}
outstruct.LCR =
{
  1x2 struct array containing the fields: % array has one element per block
    func    % name of the sbn file corresponding to each block
    section  % string parameter to be passed to the sbn files
    nextvar  % number of external variables for each element of the block
    vnmatrix % numbers of the external variables of each element
    nintvar  % number of internal variables for each element of the block
    osintvar % number of the first internal variable
    npar    % number of parameters
    nparnames% number of parameter names
    nrows    % number of rows in the block
    parnames % list of parameter names
    pvmatrix % list of parameter values for each element
}
</syntaxhighlight>}}
== File Formats ==
== File Formats ==


Line 86: Line 185:


{{Code|Model evaluator file for simple MOSFET models |<syntaxhighlight lang="octave" style="font-size:13px">
{{Code|Model evaluator file for simple MOSFET models |<syntaxhighlight lang="octave" style="font-size:13px">
function [a,b,c ,] =...
function [a, b, c] =...
func (string , m(i ,:) , extvar , intvar , t)
func (string , pvmatrix(i ,:) , extvar , intvar , t)
</syntaxhighlight>}}
</syntaxhighlight>}}
i.e. it should get as inputs:
i.e. it should get as inputs:
Line 95: Line 194:
* the current time
* the current time
and it should produce as outputs three matrices:
and it should produce as outputs three matrices:
* <math>a, b \in \mathbb{R}^{(n_extvar + n_intvar) \times (n_extvar + n_intvar)}</math>
* <math>a, b \in \mathbb{R}^{({n_{extvar}} + {n_{intvar}})  
* <math>c \in \mathbb{R}^{(n_extvar + n_intvar)}</math>
\times ({n_{extvar}} + {n_{intvar})}}</math>
* <math>c \in \mathbb{R}^{({n_{extvar}} + {n_{intvar}})}</math>
where "n_intvar" is the number of internal variables that can be assembled in the complete system matrices.
where "n_intvar" is the number of internal variables that can be assembled in the complete system matrices.


Line 473: Line 573:
</syntaxhighlight>
</syntaxhighlight>
}}
}}
=== A circuit with a linear VCVS ===
To parse an IFF format netlist of the VCVS circuit  we can use the following command
{{Code|Load the VCVS circuit structure parsing an IFF netlist |<syntaxhighlight lang="octave" style="font-size:13px">
outstruct = prs_iff ("vcvs");
</syntaxhighlight>
}}
The IFF netlist consists of the .cir file named "vcvs.cir" shown below
{{Code|IFF netlist for the VCVS circuit (.cir file)|<syntaxhighlight lang="text" style="font-size:13px">
%0.1b1
% A Simple linear VCVS example
% Input voltage sources
Mvoltagesources sinwave 2 4
1 4
Ampl      f      delay    shift
1        1      0.0      0.0
1 0
END
% VCVS
Mvcvs LIN 4 1
1 1
Gain
5.0
2 0 1 0
% Resistor
Mresistors LIN  2 1
1 1
R
1
1 2
END
</syntaxhighlight>
}}
and of the .nms file named "and.nms shown below
{{Code|IFF netlist for the VCVS circuit (.nms file)|<syntaxhighlight lang="text" style="font-size:13px">
% 0.1b1
1 V_controller
2 V_controlled
</syntaxhighlight>
}}
The implementation for the VCVS with linear gain is shown below
{{Code|Model evaluator file for simple VCVS model |<syntaxhighlight lang="octave" style="font-size:13px">
function [a,b,c] = Mvcvs (string, parameters, parameternames, extvar,
                          intvar, t)
  if isempty(intvar)
    intvar = 0;
  endif
  switch string
      ##LCR part
    case "LIN"
      for ii=1:length(parameternames)
eval([parameternames{ii} "=" num2str(parameters(ii)) ";"])
      endfor
     
      j = intvar (1);
      Vin = extvar (3) - extvar (4);
      V = Vin * Gain;
     
      a = zeros (5);     
      b = [0  0  0    0    1;
          0  0  0    0    -1;
          0  0  0    0    0;
          0  0  0    0    0;
          1 -1  -Gain  Gain  0];
     
      c = [0 0 0 0 0];
      break
    otherwise
      error (["unknown section:" string])
  endswitch
endfunction
</syntaxhighlight>
}}
[[File:VCVS_result.png|thumb| Result of the VCVS simulation]]
To run a simulation with this circuit use the following commands:
{{Code|Run a simple transient simulation with the VCVS circuit |<syntaxhighlight lang="octave" style="font-size:13px">
>> x = [0 0 0 0]';
>> t = linspace(0,1,50);
>> [out, niter] = tst_backward_euler (outstruct, x, t, 1e-6, 100, pltvars, [0 1]);
</syntaxhighlight>
}}
=== Creating a model for a memristor device ===
To demonstrate how to write a model evaluator file (SBN file), we
will discuss the simplest memristor model shown in this paper by [[User:KaKiLa| KaKiLa]] et al. (Carbajal, J. P. et al.(2015). [http://www.mitpressjournals.org/doi/abs/10.1162/NECO_a_00694?journalCode=neco#.VgFNn9tSukp Memristor models for machine learning]. Neural Computation, 27(3). Learning; Materials Science. doi:10.1162/NECO_a_00694</ref>
The device model is presented in the original paper as
<math>
\left\{
\begin{array}{l}
\dot{x} = \mu I(t)\\
H(x) = \left\{
        \begin{array}{l}
        R \qquad x \le 0\\
        R - (R - r) x \qquad 0 < x < 1\\
        r \qquad x > 1
        \end{array}       
        \right.\\
V(t) = H(x) I(t)
\end{array}
\right.
</math>
The first thing to do is to change the model from impedance to admittance form
and write the constitutive relation for the internal variable in an "implicit form"
<math>
\left\{
\begin{array}{l}
  \dfrac{1}{\mu} \dot{x} + Q(V(t), x(t)) = \dfrac{1}{\mu} \dot{x} - I(t) = 0\\
  H(x) = \left\{
        \begin{array}{l}
        R \qquad x \le 0\\
        R - (R - r) x \qquad 0 < x < 1\\
        r \qquad x > 1
        \end{array}       
        \right.\\
I(t) = \dfrac{V(t)}{H(x)}
\end{array}
\right.
</math>
It is then useful to compute the derivatives for the current and for the constitutive relation
<math>
  \dfrac{\partial I}{\partial x} = -\dfrac{H'(x)}{H(x)^2}
</math>
<math>
  \dfrac{\partial I}{\partial V} = \dfrac{1}{H}
</math>
<math>
H'(x) = \left\{
        \begin{array}{l}
        0 \qquad x \le 0\\
        - (R - r) \qquad 0 < x < 1\\
        0 \qquad x > 1
        \end{array}       
        \right.
</math>
Next we define external variables (pin voltages) for the device and coupling variables (pin currents)
  N.B. the convention in existing device models is that pin currents are
      assumed to be entering the device
<math>
\begin{array}{l}
i^{+} = I(t)\\
i^{-} = -I(t)\\
V(t) = v^{+} - v^{-}
\end{array}
</math>
Now define the local state vector
<math>
z =
\left[
\begin{array}{l}
v^{+}\\
v^{-}\\
x
\end{array}
\right]
</math>
and the local equations
<math>
a \dot{z} + c(z) = 0
</math>
and finally define the local matrices for the memristor element
<math>
a = \left[\begin{array}{c c c} 0 & 0 & 0\\ 0 & 0 & 0\\ 0 & 0 & \dfrac{1}{\mu}\end{array}\right]
</math>
<math>
c(z) = \left[\begin{array}{c}\dfrac{V}{H}\\[2mm] -\dfrac{V}{H}\end{array}\right]
</math>
<math>
b(z) = \dfrac{\partial c}{\partial z} =
\left[
\begin{array}{c c c}
\dfrac{\partial I}{\partial V} & -\dfrac{\partial I}{\partial V} & \dfrac{\partial I}{\partial x}\\[2mm]
-\dfrac{\partial I}{\partial V} & +\dfrac{\partial I}{\partial V} & -\dfrac{\partial I}{\partial x}\\[2mm]
-\dfrac{\partial I}{\partial V} & +\dfrac{\partial I}{\partial V} & -\dfrac{\partial I}{\partial x}\\[2mm]
\end{array}
\right]
</math>
The resulting model implementation is the following
{{Code| memristor model implementation|<syntaxhighlight lang="octave">
function [a,b,c] = Mmemristors (string, parameters, parameternames,
                                extvar, intvar, t)
 
  if isempty(intvar)
    intvar = 0;
  endif
  switch string       
    case "STRUKOV"
      ## NLC part
      for ii=1:length(parameternames)
eval([parameternames{ii} "=" num2str(parameters(ii)) ";"])
      endfor
     
      v1 = extvar(1);
      v2 = extvar(2);
      x  = intvar(1);
      if (x <= 0)
        H = RH;
        Hp = 0;
      elseif (x > 0 && x < 1)
        H = RH - (RH - RL) * x;
        Hp = - (RH - RL);
      else
        H = RL;
        Hp = 0;
      endif
      I  = (v1-v2) / H;
      i1 = I;
      i2 = -I;
      dIdx = - Hp / H^2;
      dIdV = 1 / H;
      a = zeros (3);
      a(3, 3) = 1/ MU;
     
      b = [dIdV  -dIdV dIdx;
          -dIdV dIdV  -dIdx;
          -dIdV dIdV  -dIdx];
     
      c = [i1 i2 i2]';
     
      break;
     
    otherwise
      error (["unknown section:" string])
  endswitch
endfunction
</syntaxhighlight>
}}
As an example let's run a simulation by applying a sinusoidal signal with varying
frequency
{{Code| memristor model implementation|<syntaxhighlight lang="octave">
t = linspace (0, 1, 100);
y = [(sin (2 * pi * t * 1.5)), (sin (2 * pi * t * 4))(2:end)];
t = [t, (linspace (1, 2, 100))(2:end)];
pwl = [t; y](:).';
c1.LCR = [];
c1.totextvar = 1;
c1.NLC(1).("func") = "Mvoltagesources";
c1.NLC(1).("section") = "pwl";
c1.NLC(1).("nextvar") = 2;
c1.NLC(1).("npar") = 398;
c1.NLC(1).("nrows") = 1;
c1.NLC(1).("nparnames") = 0;
c1.NLC(1).("parnames") = {};
c1.NLC(1).("pvmatrix") = pwl;
c1.NLC(1).("vnmatrix") = [1 0];
c1.NLC(1).("nintvar") = 1;
c1.NLC(1).("osintvar") = 0;
c1.NLC(2).("func") = "Mmemristors";
c1.NLC(2).("section") =  "STRUKOV";
c1.NLC(2).("nextvar") =  2;
c1.NLC(2).("npar") =  3;
c1.totintvar = 2;
c1.namesn = [1  2  3];
c1.NLC(2).("nrows") =  1;
c1.NLC(2).("nparnames") =  3;
c1.NLC(2).("parnames") =  {"MU", "RH", "RL"};
c1.NLC(2).("pvmatrix") =  [1.8e3 1e3 1];
c1.NLC(2).("vnmatrix") =  [1 0];
c1.NLC(2).("nintvar") =  1;
c1.NLC(2).("osintvar") =  1;
c1.namess = {"voltage", "current", "x"};
start = [0 0 0.1];
%% c = prs_iff ("memristor_example");
out = tst_backward_euler (c1, start.', t, 1e-2, 300, {'voltage', 'current'});
plot (out(1, 1:100), out (2, 1:100), 'r',
      out(1, 101:end), out(2, 101:end), 'k')
legend ("frequency 1 Hz", "frequency 4 Hz")
</syntaxhighlight>
}}
[[File:memristor.png|thumb| Memristor simulation result]]
The results are shown in the figure to the right.
[[Category:Octave Forge]]
Retrieved from "https://wiki.octave.org/Special:MobileDiff/6637...11181"

Navigation menu