1 ## Copyright (C) 2012 Marco Merlin
3 ## This file is part of:
4 ## OCS - A Circuit Simulator for Octave
6 ## OCS is free software; you can redistribute it and/or modify
7 ## it under the terms of the GNU General Public License as published by
8 ## the Free Software Foundation; either version 2 of the License, or
9 ## (at your option) any later version.
11 ## OCS is distributed in the hope that it will be useful,
12 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
13 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 ## GNU General Public License for more details.
16 ## You should have received a copy of the GNU General Public License
17 ## along with OCS; If not, see <http://www.gnu.org/licenses/>.
19 ## author: Marco Merlin <marcomerli _AT_ gmail.com>
20 ## based on prs_iff which is (C) Carlo de Falco and Massimiliano Culpo
23 ## @deftypefn {Function File} {[@var{stuct}]} = prs_spice (@var{filename})
25 ## Circuit file parser that can interpret a subset of the spice file format.
27 ## @code{prs_spice} currently supports the following set of "Element Cards"
36 ## Cname anode knode modelname <parameters>
41 ## Mname gnode dnode snode bnode modelname <parameters>
44 ## N.B.: one instance of a MOS element MUST be preceeded (everywhere in the file) by the declaration of the related model.
47 ## .MODEL mynmos NMOS( k=1e-4 Vth=0.1 rd=1e6)
48 ## M2 Vgate 0 Vdrain 0 mynmos
56 ## @item Voltage sources:
58 ## Vname n+ n- <dcvalue> <transvalue>
61 ## Transvalue specifies a transient voltage source
63 ## SIN(VO VA FREQ TD THETA)
68 ## @item VA (amplitude)
69 ## @item FREQ (frequency)
71 ## @item THETA (damping factor)
77 ## VO + VA*exp(-(time-TD)*THETA)*sine(twopi*FREQ*(time+TD))
80 ## Currently the damping factor has no effect.
84 ## PULSE(V1 V2 TD TR TF PW PER)
89 ## @item V1 (initial value)
90 ## @item V2 (pulsed value)
91 ## @item TD (delay time)
92 ## @item TR (rise time)
93 ## @item TF (fall time)
94 ## @item PW (pulse width)
98 ## Currently rise and fall time are not implemented yet.
100 ## @item .MODEL cards
101 ## Defines a model for semiconductor devices
104 ## .MODEL MNAME TYPE(PNAME1=PVAL1 PNAME2=PVAL2 ... )
109 ## @item NMOS N-channel MOSFET model
110 ## @item PMOS P-channel MOSFET model
111 ## @item D diode model
114 ## The parameter "LEVEL" is currently assigned to the field "section" in the call
115 ## of the element functions by the solver.
116 ## Currently supported values for the parameter LEVEL for NMOS and PMOS are:
121 ## (see documentation of function Mdiode).
123 ## Currently supported values for the parameter LEVEL for D are:
127 ## (see documentation of functions Mnmosfet and Mpmosfet).
130 ## @seealso{prs_iff,Mdiode,Mnmosfet,Mpmosfet}
133 function outstruct = prs_spice (name)
136 if (nargin != 1 || !ischar (name))
137 error ("prs_spice: wrong input.")
142 outstruct = struct ("LCR", [],
146 global ndsvec; # Vector of circuit nodes
152 outstruct.totintvar = 0; # Number of internal variables
154 count = struct ("NLC", struct ("n", 0,
156 "LCR", struct ("n", 0,
160 models_list = struct ("mname", {""},
167 filename = [name ".spc"];
168 if (exist (filename) != 2) ## isempty (file_in_path (".", filename))
169 error (["prs_spice: .spc file not found:" filename]);
171 fid = fopen (filename, "r");
178 line = strtrim (line);
180 %% exclude empty lines
183 if (! strncmpi (line, '*', 1))
184 %% lines here aren't comments
185 if (strncmpi (line, '+', 1))
186 %% this line has to be concatenated to previous one
188 fullline = [fullline line];
190 %% these lines are not a concatenation
192 %% line echo for debug
195 %% compute fullline here!
196 %#[outstruct, intvar, count] = lineParse (upper (fullline), outstruct, count, intvar);
198 %# NB: case-sensitive because of parameter names
199 %#[outstruct, intvar, count] = lineParse (fullline, outstruct, intvar, count);
200 [outstruct, count] = lineParse (fullline, outstruct, count);
203 end %if (strncmpi (line, '+', 1))
204 end %if (~strncmpi (line, '*', 1))
205 end % if length (line)
207 lineCounter = lineCounter+1;
208 end % while ~feof (fid)
212 ## NB: case-sensitive because of parameter names
213 %#[outstruct, intvar, count] = lineParse (upper (fullline), outstruct, count, intvar);
214 %#[outstruct, intvar, count] = lineParse (fullline, outstruct, intvar, count);
215 [outstruct, count] = lineParse (fullline, outstruct, count);
219 error ('Input file not found!');
223 ## Set the number of internal and external variables
224 nnodes = length (unique (ndsvec));
225 maxidx = max (ndsvec);
227 if (nnodes <= (maxidx+1))
228 ## If the valid file is a subcircuit it may happen
229 ## that nnodes == maxidx, otherwise nnodes == (maxidx+1)
230 outstruct.totextvar = max (ndsvec);
232 error ("prs_spice: hanging nodes in circuit %s", name);
235 ## set node names as variable names
236 for ii = 1:length (nodes_list)
237 outstruct.namesn (ii) = ii-1;
239 outstruct.namess = horzcat (nodes_list, intvar_list);
242 ##outstruct.totintvar
245 ## NLC block intvar count update
246 function outstruct = NLCintvar (outstruct, NLCcount, name)
251 ## set node names for NLC subcircuit
252 ##for NLCcount = 1:count.NLC.n;
253 ##NLCcount = count.NLC.n;
254 ## set vnmatrix for NLC subcircuit
256 outstruct.NLC(NLCcount).vnmatrix(:)];
257 ## Compute internal variables cycling over each
258 ## element in the section
259 ##for iel = 1:outstruct.NLC(NLCcount).nrows
260 iel = outstruct.NLC(NLCcount).nrows;
261 [a, b, c] = feval (outstruct.NLC(NLCcount).func,
262 outstruct.NLC(NLCcount).section,
263 outstruct.NLC(NLCcount).pvmatrix(iel, :),
264 outstruct.NLC(NLCcount).parnames,
265 zeros (outstruct.NLC(NLCcount).nextvar, 1),
269 ## FIXME: if all the element in the same section share the
270 ## same number of internal variables, the for cycle can be
271 ## substituted by only one call
272 outstruct.NLC(NLCcount).nintvar(iel) = columns (a) - outstruct.NLC(NLCcount).nextvar;
273 outstruct.NLC(NLCcount).osintvar(iel) = outstruct.totintvar;
275 outstruct.totintvar += outstruct.NLC(NLCcount).nintvar(iel);
276 if outstruct.NLC(NLCcount).nintvar(iel)>0
277 intvar_list{outstruct.totintvar} = ["I(" name ")"];
280 ##endfor # NLCcount = 1:count.NLC.n;
283 ## LCR block intvar count update
284 function outstruct = LCRintvar (outstruct, LCRcount, name)
287 ## set node names for LCR subcircuit
288 ## for LCRcount = 1:count.LCR.n;
289 ## LCRcount = count.LCR.n;
290 ## set vnmatrix for LCR subcircuit
292 outstruct.LCR(LCRcount).vnmatrix(:)];
294 ## Compute internal variables cycling over each
295 ## element in the section
296 ##for iel = 1:outstruct.LCR(LCRcount).nrows
297 iel = outstruct.LCR(LCRcount).nrows;
298 [a, b, c] = feval (outstruct.LCR(LCRcount).func,
299 outstruct.LCR(LCRcount).section,
300 outstruct.LCR(LCRcount).pvmatrix(iel, :),
301 outstruct.LCR(LCRcount).parnames,
302 zeros(outstruct.LCR(LCRcount).nextvar, 1),
306 ## FIXME: if all the element in the same section share the
307 ## same number of internal variables, the for cycle can be
308 ## substituted by only one call
309 outstruct.LCR(LCRcount).nintvar(iel) = columns (a) - outstruct.LCR(LCRcount).nextvar;
310 ##outstruct.LCR(LCRcount).osintvar(iel) = intvar;
311 outstruct.LCR(LCRcount).osintvar(iel) = outstruct.totintvar;
313 ##intvar += outstruct.LCR(LCRcount).nintvar(iel);
314 outstruct.totintvar += outstruct.LCR(LCRcount).nintvar(iel);
315 if outstruct.LCR(LCRcount).nintvar(iel)>0
316 intvar_list{outstruct.totintvar} = ["I(" name ")"];
319 ## endfor # LCRcount = 1:count.LCR.n;
322 ## Parses a single line
323 function [outstruct, count] = lineParse (line, outstruct, count)
328 [outstruct, count] = prs_spice_C (line, outstruct, count);
330 [outstruct, count] = prs_spice_D (line, outstruct, count);
332 ##mspINP2E (line, lineCounter);
334 ##mspINP2F (line, lineCounter);
336 ##mspINP2G (line, lineCounter);
338 ##mspINP2H (line, lineCounter);
340 ##mspINP2I (line, lineCounter);
343 ##mspINP2K (line, lineCounter);
345 ##mspINP2L (line, lineCounter);
347 ## FIXME: just for nMOS devices!
348 [outstruct, count] = prs_spice_M (line, outstruct, count);
350 ## temporarily assigned to pMOS devices.
351 ##[outstruct, count] = prs_spice_P (line, outstruct, count);
354 [outstruct, count] = prs_spice_R (line, outstruct, count);
359 [outstruct, count] = prs_spice_V (line, outstruct, count);
364 [outstruct, count] = prs_spice_dot (line, outstruct, count);
366 warn = sprintf (['prs_spice: Unsupported circuit element in line: ' line ]);
368 end % switch (line (1))
369 end %if length (line)
373 ## adds an NLC element to outstruct
374 function [outstruct, count] = addNLCelement (outstruct, count, element, section, nextvar, npar, nparnames, parnames, vnmatrix, pvmatrix)
375 ## check whether the element type already exists in output structure
376 ## the search key is (func, section)
377 ##[tf, idx] = ismember ({element}, count.NLC.list);
379 ## the element still does not exist
382 count.NLC.list{count.NLC.n} = element;
384 ## if the element doesn't exists, add it to the output structure
385 outstruct.NLC(count.NLC.n).func = element;
386 outstruct.NLC(count.NLC.n).section = section;
387 outstruct.NLC(count.NLC.n).nextvar = nextvar;
388 outstruct.NLC(count.NLC.n).npar = npar;
389 outstruct.NLC(count.NLC.n).nrows = 1;
390 outstruct.NLC(count.NLC.n).nparnames = nparnames;
391 outstruct.NLC(count.NLC.n).parnames = parnames;
392 outstruct.NLC(count.NLC.n).vnmatrix = vnmatrix;
393 outstruct.NLC(count.NLC.n).pvmatrix = pvmatrix;
396 ## for ii = 1:length (idx)
397 ## if strcmp (outstruct.NLC(idx(ii)).section, section)
400 ## endif #!strcmp (outstruct.NLC(idx(ii)).section, section)
401 ## endfor #ii = 1:length (idx)
404 ## the section does not exist
406 ## the element still does not exist
409 ## count.NLC.list{count.NLC.n} = element;
411 ## if the element doesn't exists, add it to the output structure
412 ## outstruct.NLC(count.NLC.n).func = element;
413 ## outstruct.NLC(count.NLC.n).section = section;
414 ## outstruct.NLC(count.NLC.n).nextvar = nextvar;
415 ## outstruct.NLC(count.NLC.n).npar = npar;
416 ## outstruct.NLC(count.NLC.n).nrows = 1;
417 ## outstruct.NLC(count.NLC.n).nparnames = nparnames;
418 ## outstruct.NLC(count.NLC.n).parnames = parnames;
419 ## outstruct.NLC(count.NLC.n).vnmatrix = vnmatrix;
420 ## outstruct.NLC(count.NLC.n).pvmatrix = pvmatrix;
423 ## the couple (element, section) already exists, so add a row in the structure
424 ## add an element to the structure
425 ## outstruct.NLC(idx(ii)).nrows++;
427 ## update parameter value and connectivity matrix
428 ## [outstruct.NLC(idx(ii)).vnmatrix] = [outstruct.NLC(idx(ii)).vnmatrix; vnmatrix];
429 ## [outstruct.NLC(idx(ii)).pvmatrix] = [outstruct.NLC(idx(ii)).pvmatrix; pvmatrix];
434 ## adds an LCR element to outstruct
435 function [outstruct, count] = addLCRelement (outstruct, count, element, section, nextvar, npar, nparnames, parnames, vnmatrix, pvmatrix)
436 ## check whether the element type already exists in output structure
437 ## the search key is (func, section)
438 [tf, idx] = ismember ({element}, count.LCR.list);
440 ## the element still does not exist
443 count.LCR.list{count.LCR.n} = element;
445 ## if the element doesn't exists, add it to the output structure
446 outstruct.LCR(count.LCR.n).func = element;
447 outstruct.LCR(count.LCR.n).section = section;
448 outstruct.LCR(count.LCR.n).nextvar = nextvar;
449 outstruct.LCR(count.LCR.n).npar = npar;
450 outstruct.LCR(count.LCR.n).nrows = 1;
451 outstruct.LCR(count.LCR.n).nparnames = nparnames;
452 outstruct.LCR(count.LCR.n).parnames = parnames;
453 outstruct.LCR(count.LCR.n).vnmatrix = vnmatrix;
454 outstruct.LCR(count.LCR.n).pvmatrix = pvmatrix;
457 for ii = 1:length (idx)
458 if strcmp (outstruct.LCR(idx(ii)).section, section)
461 endif #!strcmp (outstruct.LCR(idx(ii)).section, section)
462 endfor #ii = 1:length (idx)
465 ## the section does not exist
467 ## the element still does not exist
470 count.LCR.list{count.LCR.n} = element;
472 ## if the element doesn't exists, add it to the output structure
473 outstruct.LCR(count.LCR.n).func = element;
474 outstruct.LCR(count.LCR.n).section = section;
475 outstruct.LCR(count.LCR.n).nextvar = nextvar;
476 outstruct.LCR(count.LCR.n).npar = npar;
477 outstruct.LCR(count.LCR.n).nrows = 1;
478 outstruct.LCR(count.LCR.n).nparnames = nparnames;
479 outstruct.LCR(count.LCR.n).parnames = parnames;
480 outstruct.LCR(count.LCR.n).vnmatrix = vnmatrix;
481 outstruct.LCR(count.LCR.n).pvmatrix = pvmatrix;
484 ## the couple (element, section) already exists, so add a row in the structure
485 ## add an element to the structure
486 outstruct.LCR(idx(ii)).nrows++;
488 ## update parameter value and connectivity matrix
489 [outstruct.LCR(idx(ii)).vnmatrix] = [outstruct.LCR(idx(ii)).vnmatrix; vnmatrix];
490 [outstruct.LCR(idx(ii)).pvmatrix] = [outstruct.LCR(idx(ii)).pvmatrix; pvmatrix];
495 ## converts a blank separated values string into a cell array
496 function ca = str2ca (str)
497 ca = regexpi (str, '[^ \s=]+', 'match');
500 ## replaces the tokens c_old with c_new in string inString
501 function outString = strReplace (inString, c_old, c_new)
502 outString = inString;
504 for idx = 1:length (c_old)
506 outString = strrep (outString, c_old{idx}, c_new{idx});
511 ## returns the numeric value of a string
512 function num = literal2num (string)
513 literals = {'MEG' 'MIL' 'A' 'F' 'P' 'N' 'U' 'M' 'K' 'G' 'T'};
514 numerics = {'e6' '*25.4e-6' 'e-18' 'e-15' 'e-12' 'e-9' 'e-6' 'e-3' 'e3' 'e9' 'e12'};
515 newstr = strReplace (upper (string), literals, numerics);
516 num = str2num (newstr);
519 function syntaxError (line)
520 warnstr = sprintf ("Syntax error in line: %s", line);
524 function [outstruct, count] = prs_spice_C (line, outstruct, count)
525 element = "Mcapacitors";
526 ## check wheter the element type already exists in output structure
527 [tf, idx] = ismember ({element}, count.NLC.list);
531 count.NLC.list{count.NLC.n} = element;
534 ## if the element doesn't exists, add it to the output structure
535 outstruct.NLC(count.NLC.n).func = element;
536 outstruct.NLC(count.NLC.n).section = "LIN";
537 outstruct.NLC(count.NLC.n).nextvar = 2;
538 outstruct.NLC(count.NLC.n).npar = 1;
539 outstruct.NLC(count.NLC.n).nrows = 0;
540 outstruct.NLC(count.NLC.n).nparnames = 1;
541 outstruct.NLC(count.NLC.n).parnames = {"C"};
542 outstruct.NLC(count.NLC.n).vnmatrix = [];
543 outstruct.NLC(count.NLC.n).pvmatrix = [];
546 ## add an element to the structure
547 ##outstruct.NLC(idx).nrows++;
548 outstruct.NLC(count.NLC.n).nrows++;
550 ## convert input line string into cell array
554 ## update parameter value and connectivity matrix
555 ##[outstruct.NLC(idx).vnmatrix] = [outstruct.NLC(idx).vnmatrix; add_nodes(ca(2:3))];
556 ##[outstruct.NLC(idx).pvmatrix] = [outstruct.NLC(idx).pvmatrix; literal2num(ca{4})];
557 [outstruct.NLC(count.NLC.n).vnmatrix] = [outstruct.NLC(count.NLC.n).vnmatrix; add_nodes(ca(2:3))];
558 [outstruct.NLC(count.NLC.n).pvmatrix] = [outstruct.NLC(count.NLC.n).pvmatrix; literal2num(ca{4})];
563 ##outstruct = NLCintvar (outstruct, idx, ca{1});
564 outstruct = NLCintvar(outstruct, count.NLC.n, ca{1});
567 function [outstruct, count] = prs_spice_D (line, outstruct, count)
569 ## check wheter the element type already exists in output structure
570 ##[tf, idx] = ismember ({element}, count.NLC.list);
574 count.NLC.list{count.NLC.n} = element;
577 ## if the element doesn't exists, add it to the output structure
578 outstruct.NLC(count.NLC.n).func = element;
579 outstruct.NLC(count.NLC.n).section = "simple";
580 outstruct.NLC(count.NLC.n).nextvar = 2;
581 outstruct.NLC(count.NLC.n).npar = 0;
582 outstruct.NLC(count.NLC.n).nrows = 1;
583 outstruct.NLC(count.NLC.n).nparnames = 0;
584 outstruct.NLC(count.NLC.n).parnames = {};
585 outstruct.NLC(count.NLC.n).vnmatrix = [];
586 outstruct.NLC(count.NLC.n).pvmatrix = [];
589 ## convert input line string into cell array
591 ## update parameter value and connectivity matrix
592 ##[outstruct.NLC(idx).vnmatrix] = [outstruct.NLC(idx).vnmatrix; add_nodes(ca(2:3))];
593 [outstruct.NLC(count.NLC.n).vnmatrix] = [outstruct.NLC(count.NLC.n).vnmatrix; add_nodes(ca(2:3))];
595 ## check if parameters are specified
596 for prm_idx = 1:(length (ca)-3)/2
597 ## [tf, str_idx] = ismember (outstruct.NLC(count.NLC.n).parnames{prm_idx}, ca);
598 ## TODO: can the loop be executed in a single operation?
599 ## [tf, str_idx] = ismember (outstruct.NLC(count.NLC.n).parnames, ca);
600 if length (ca) >= 1+prm_idx*2
601 ## find specified parameter name
604 outstruct.NLC(count.NLC.n).npar++;
605 outstruct.NLC(count.NLC.n).nparnames++;
606 outstruct.NLC(count.NLC.n).parnames{1, prm_idx} = ca{2*prm_idx+2};
607 outstruct.NLC(count.NLC.n).pvmatrix(1, prm_idx) = literal2num (ca{2*prm_idx+3});
609 ## TODO: set to a default value undefined parameters, instead of rising an error?
610 ## errstr = sprintf ("Undefined parameter %s in line: %s", outstruct.NLC(count.NLC.n).parnames{prm_idx}, line);
618 outstruct = NLCintvar(outstruct, count.NLC.n, ca{1});
621 function [outstruct, count] = prs_spice_M (line, outstruct, count)
623 ##element = "Mnmosfet";
624 ## check wheter the element type already exists in output structure
625 ##[tf, idx] = ismember ({element}, count.NLC.list);
629 ##count.NLC.list{count.NLC.n} = element;
632 ## if the element doesn't exists, add it to the output structure
633 ##outstruct.NLC(count.NLC.n).func = element;
634 ##outstruct.NLC(count.NLC.n).section = "simple";
635 ## outstruct.NLC(count.NLC.n).section = "lincap";
636 outstruct.NLC(count.NLC.n).nextvar = 4;
637 outstruct.NLC(count.NLC.n).npar = 0;
638 outstruct.NLC(count.NLC.n).nrows = 1;
639 outstruct.NLC(count.NLC.n).nparnames = 0;
640 outstruct.NLC(count.NLC.n).parnames = {};
641 ##outstruct.NLC(count.NLC.n).vnmatrix = [];
642 outstruct.NLC(count.NLC.n).pvmatrix = [];
645 ## convert input line string into cell array
647 ## update parameter value and connectivity matrix
648 ##[outstruct.NLC(idx).vnmatrix] = [add_nodes(ca(2:5))];
649 [outstruct.NLC(count.NLC.n).vnmatrix] = [add_nodes(ca(2:5))];
651 [tf, idx] = ismember (ca{6}, models_list.mname);
654 outstruct.NLC(count.NLC.n).func = models_list.melement{idx};
655 outstruct.NLC(count.NLC.n).section = models_list.msection{idx};
657 ## check if parameters are specified
658 for prm_idx = 1:(length (ca)-6)/2
659 if length (ca)>=4+prm_idx*2
660 outstruct.NLC(count.NLC.n).npar++;
661 outstruct.NLC(count.NLC.n).nparnames++;
662 outstruct.NLC(count.NLC.n).parnames{1, prm_idx} = ca{2*prm_idx+5};
663 outstruct.NLC(count.NLC.n).pvmatrix(1, prm_idx) = literal2num (ca{2*prm_idx+6});
669 ## add model parameters to list
670 prm_idx = outstruct.NLC(count.NLC.n).npar;
671 len = length (models_list.plist{idx}.pnames);
673 outstruct.NLC(count.NLC.n).parnames{prm_idx+mpidx} = models_list.plist{idx}.pnames{mpidx};
674 outstruct.NLC(count.NLC.n).pvmatrix(1, prm_idx+mpidx) = models_list.plist{idx}.pvalues(mpidx);
676 outstruct.NLC(count.NLC.n).npar = outstruct.NLC(count.NLC.n).npar+len;
677 outstruct.NLC(count.NLC.n).nparnanames = outstruct.NLC(count.NLC.n).nparnames+len;
679 outstruct = NLCintvar(outstruct, count.NLC.n, ca{1});
685 function [outstruct, count]= prs_spice_P (line, outstruct, count)
686 element = "Mpmosfet";
687 ## check wheter the element type already exists in output structure
690 count.NLC.list{count.NLC.n} = element;
692 ## if the element doesn't exists, add it to the output structure
693 outstruct.NLC(count.NLC.n).func = element;
694 outstruct.NLC(count.NLC.n).section = "simple";
695 ## outstruct.NLC(count.NLC.n).section = "lincap";
696 outstruct.NLC(count.NLC.n).nextvar = 4;
697 outstruct.NLC(count.NLC.n).npar = 0;
698 outstruct.NLC(count.NLC.n).nrows = 1;
699 outstruct.NLC(count.NLC.n).nparnames = 0;
700 outstruct.NLC(count.NLC.n).parnames = {};
701 ##outstruct.NLC(count.NLC.n).vnmatrix = [];
702 outstruct.NLC(count.NLC.n).pvmatrix = [];
704 ## convert input line string into cell array
706 ## update parameter value and connectivity matrix
707 ##[outstruct.NLC(idx).vnmatrix] = [add_nodes(ca(2:5))];
708 [outstruct.NLC(count.NLC.n).vnmatrix] = [add_nodes(ca(2:5))];
710 ## check if parameters are specified
711 for prm_idx = 1:(length (ca)-5)/2
712 ## [tf, str_idx] = ismember (outstruct.NLC(count.NLC.n).parnames{prm_idx}, ca);
713 ## TODO: can the loop be executed in a single operation?
714 ## [tf, str_idx] = ismember (outstruct.NLC(count.NLC.n).parnames, ca);
715 if (length (ca) >= 3+prm_idx*2)
716 ## find specified parameter name
719 outstruct.NLC(count.NLC.n).npar++;
720 outstruct.NLC(count.NLC.n).nparnames++;
721 outstruct.NLC(count.NLC.n).parnames{1, prm_idx} = ca{2*prm_idx+4};
722 outstruct.NLC(count.NLC.n).pvmatrix(1, prm_idx) = literal2num (ca{2*prm_idx+5});
724 ## TODO: set to a default value undefined parameters, instead of rising an error?
725 ## errstr=sprintf("Undefined parameter %s in line: %s", outstruct.NLC(count.NLC.n).parnames{prm_idx}, line);
733 outstruct = NLCintvar(outstruct, count.NLC.n, ca{1});
737 function [outstruct, count]= prs_spice_R (line, outstruct, count)
739 element = "Mresistors";
740 ## check wheter the element type already exists in output structure
741 [tf, idx] = ismember ({element}, count.LCR.list);
745 count.LCR.list{count.LCR.n} = element;
748 ## if the element doesn't exists, add it to the output structure
749 outstruct.LCR(count.LCR.n).func = element;
750 outstruct.LCR(count.LCR.n).section = "LIN";
751 outstruct.LCR(count.LCR.n).nextvar = 2;
752 outstruct.LCR(count.LCR.n).npar = 1;
753 outstruct.LCR(count.LCR.n).nrows = 0;
754 outstruct.LCR(count.LCR.n).nparnames = 1;
755 outstruct.LCR(count.LCR.n).parnames = {"R"};
756 outstruct.LCR(count.LCR.n).vnmatrix = [];
757 outstruct.LCR(count.LCR.n).pvmatrix = [];
760 ## add an element to the structure
761 ##outstruct.LCR(idx).nrows++;
762 outstruct.LCR(count.LCR.n).nrows++;
764 ## convert input line string into cell array
768 ## update parameter value and connectivity matrix
769 ##[outstruct.LCR(idx).vnmatrix] = [outstruct.LCR(idx).vnmatrix; add_nodes(ca(2:3))];
770 ##[outstruct.LCR(idx).pvmatrix] = [outstruct.LCR(idx).pvmatrix; literal2num(ca{4})];
771 [outstruct.LCR(count.LCR.n).vnmatrix] = [outstruct.LCR(count.LCR.n).vnmatrix; add_nodes(ca(2:3))];
772 [outstruct.LCR(count.LCR.n).pvmatrix] = [outstruct.LCR(count.LCR.n).pvmatrix; literal2num(ca{4})];
777 outstruct = LCRintvar(outstruct, count.LCR.n, ca{1});
780 function [outstruct, count] = prs_spice_V (line, outstruct, count)
783 ## SIN(VO VA FREQ TD THETA)
789 ## THETA (damping factor)
792 ## TD to TSTOP: VO + VA*exp(-(time-TD)*THETA)*sine(twopi*FREQ*(time+TD)
794 ## check if it is a sinwave generator
795 sine = regexp (line, '(?<stim>SIN)[\s]*\((?<prms>.+)\)', 'names');
798 ## PULSE(V1 V2 TD TR TF PW PER)
800 ## parameters default values units
801 ## V1 (initial value) Volts or Amps
802 ## V2 (pulsed value) Volts or Amps
803 ## TD (delay time) 0.0 seconds
804 ## TR (rise time) TSTEP seconds
805 ## TF (fall time) TSTEP seconds
806 ## PW (pulse width) TSTOP seconds
807 ## PER (period) TSTOP seconds
809 ## check if it is a pulse generator
810 pulse = regexp (line, '(?<stim>PULSE)[\s]*\((?<prms>.+)\)', 'names');
812 if (! isempty (sine.stim))
814 ca = str2ca (sine.prms);
815 if (length (ca) == 5)
816 vo = literal2num (ca{1});
817 va = literal2num (ca{2});
818 freq = literal2num (ca{3});
819 td = literal2num (ca{4});
820 theta = literal2num (ca{5});
822 pvmatrix = [va freq td vo];
824 element = "Mvoltagesources";
830 parnames = {"Ampl", "f", "delay", "shift"};
831 ## convert input line string into cell array
833 vnmatrix = add_nodes (ca(2:3));
834 [outstruct, count] = addNLCelement (outstruct, count, element, section, nextvar, npar, nparnames, parnames, vnmatrix, pvmatrix);
835 outstruct = NLCintvar(outstruct, count.NLC.n, ca{1});
836 else #length(ca) == 5
838 endif #length (ca) == 5
839 elseif (! isempty (pulse.stim))
840 ca = str2ca (pulse.prms);
842 low = literal2num (ca{1});
843 high = literal2num (ca{2});
844 delay = literal2num (ca{3});
845 ## TODO: add rise and fall times!
846 ## tr = literal2num (ca{4});
847 ## tf = literal2num (ca{5});
848 thigh = literal2num (ca{6});
849 period = literal2num (ca{7})
852 pvmatrix = [low high tlow thigh delay];
854 element = "Mvoltagesources";
855 section = "squarewave";
860 parnames = {"low", "high", "tlow", "thigh", "delay"};
861 ## convert input line string into cell array
863 vnmatrix = add_nodes (ca(2:3));
864 [outstruct, count] = addNLCelement (outstruct, count, element, section, nextvar, npar, nparnames, parnames, vnmatrix, pvmatrix);
865 outstruct = NLCintvar(outstruct, count.NLC.n, ca{1});
866 else ##length (ca) == 7
868 endif ##length (ca) == 7
869 else ##~isempty (tran.stim)
871 element = "Mvoltagesources";
878 ## convert input line string into cell array
880 vnmatrix = add_nodes (ca(2:3));
882 pvmatrix = literal2num (ca{4});
883 [outstruct, count] = addLCRelement (outstruct, count, element, section, nextvar, npar, nparnames, parnames, vnmatrix, pvmatrix);
884 outstruct = LCRintvar(outstruct, count.LCR.n, ca{1});
885 endif #~isempty(tran.stim)
888 function [outstruct, count] = prs_spice_dot (line, outstruct, count)
889 ## .MODEL MNAME TYPE(PNAMEl = PVALl PNAME2 = PVAL2 ... )
894 ## URC Uniform Distributed RC model
898 ## NJF N-channel JFET model
899 ## PJF P-channel lFET model
900 ## NMOS N-channel MOSFET model
901 ## PMOS P-channel MOSFET model
902 ## NMF N-channel MESFET model
903 ## PMF P-channel MESFET model
904 ## SW voltage controlled switch
905 ## CSW current controlled switch
909 model = regexp (line, '.MODEL[\s]+(?<mname>[\S]+)[\s]+(?<mtype>R|C|URC|D|NPN|PNP|NJF|PJF|NMOS|PMOS|NMF|PMF|SW|CSW)[\s]*\([\s]*(?<prms>.+)[\s]*\)', 'names');
922 element = "Mnmosfet";
924 element = "Mpmosfet";
933 ## get model level (=section)
934 level = regexp (model.prms, 'LEVEL=(?<level>[\S]+)[\s]+(?<prms>.+)', 'names');
935 if isempty (level.level)
938 section = level.level;
939 model.prms = level.prms;
942 midx = length (models_list.mname)+1;
944 models_list.mname{midx} = model.mname;
945 models_list.melement{midx} = element;
946 models_list.msection{midx} = section;
948 ## set model parameters
949 ca = str2ca (model.prms);
950 midx = length (models_list.mname);
952 models_list.plist{midx} = struct ("pnames", {""}, "pvalues", []);
954 for prm_idx = 1:length (ca)/2
955 ## save parameter name
956 models_list.plist{midx}.pnames{prm_idx} = ca{2*prm_idx-1};
957 ## save parameter value
958 models_list.plist{midx}.pvalues = [models_list.plist{midx}.pvalues literal2num(ca{2*prm_idx})];
961 endif #!isempty (model)
964 ## adds nodes to nodes_list and returns the indexes of the added nodes
965 function indexes = add_nodes (nodes)
967 for ii = 1:length (nodes)
968 [tf, idx] = ismember (nodes(ii), nodes_list);
972 indexes(ii) = length (nodes_list);
973 nodes_list{length (nodes_list)+1} = nodes{ii};
979 %! outstruct = prs_spice ("rlc");
980 %! vin = linspace (0, 10, 10);
981 %! x = zeros (outstruct.totextvar+outstruct.totintvar, 1);
983 %! for idx = 1:length (vin)
984 %! outstruct.NLC(1).pvmatrix(1) = vin(idx);
985 %! out = nls_stationary (outstruct, x, 1e-15, 100);
986 %! vout(idx) = out(2);
995 %! cir = menu ("Chose the circuit to analyze:",
996 %! "AND (Simple algebraic MOS-FET model)",
997 %! "AND (Simple MOS-FET model with parasitic capacitances)",
998 %! "Diode clamper (Simple exponential diode model)",
999 %! "CMOS-inverter circuit (Simple algebraic MOS-FET model)",
1000 %! "n-MOS analog amplifier (Simple algebraic MOS-FET model)",
1001 %! "Linear RC circuit",
1002 %! "Diode bridge rectifier",
1007 %! outstruct = prs_spice ("and");
1008 %! x = [.5 .5 .33 .66 .5 1 0 0 1 ]';
1009 %! t = linspace (0, .5, 100);
1010 %! pltvars = {"Va", "Vb", "Va_and_b"};
1015 %! outstruct = prs_spice ("and2");
1016 %! x = [.8;.8;0.00232;0.00465;.8;
1017 %! .8;0.00232;0.00465;0.00000;
1018 %! 0.0;0.0;0.0;0.00232;0.0;
1019 %! 0.0;0.0;0.0;1.0;0.0;-0.0;
1020 %! 0.0;1.0;0.00465;0.0;0.0;
1021 %! -0.0;1.0;0.00465;0.0;
1022 %! 0.0;0.0;1.0;1.0;0.0;0.0;0.0;
1024 %! t = linspace (.25e-6, .5e-6, 100);
1026 %! pltvars = {"Va", "Vb", "Va_and_b"};
1030 %! outstruct = prs_spice ("diode");
1031 %! x = [0 0 0 0 0]';
1032 %! t = linspace (0, 3e-10, 200);
1034 %! pltvars = {"Vin", "Vout", "V2"};
1038 %! outstruct = prs_spice ("inverter");
1039 %! x = [.5 .5 1 0 0]';
1040 %! t = linspace (0, 1, 100);
1042 %! pltvars={"Vgate", "Vdrain"};
1046 %! outstruct = prs_spice ("nmos");
1047 %! x = [1 .03 1 0 0]';
1048 %! t = linspace (0, 1, 50);
1050 %! pltvars = {"Vgate", "Vdrain"};
1054 %! outstruct = prs_spice ("rcs");
1056 %! t = linspace (0, 2e-5, 100);
1058 %! pltvars = {"Vout"};
1062 %! outstruct = prs_spice ("rect");
1064 %! t = linspace (0, 3e-10, 60);
1066 %! pltvars = {"Vin", "Voutlow", "Vouthigh"};
1070 %! outstruct = prs_spice ("rlc")
1071 %! #x = [0 0 0 0 0]';
1074 %! x = [0 0 0 0 0 0 ]';
1075 %! t = linspace (0, 2e-5, 200);
1077 %! #pltvars = {"Vin", "Vout"};
1078 %! pltvars = {"I(C3)"};
1079 %! #pltvars = {"Vout"};
1083 %! error ("unknown circuit");
1087 %! slv = menu("Chose the solver to use:",
1097 %! out = zeros (rows (x), columns (t));
1101 %! out = tst_backward_euler (outstruct, x, t, tol, maxit, pltvars);
1102 %! # out = TSTbweuler (outstruct, x, t, tol, maxit, pltvars);
1104 %! out = tst_daspk (outstruct, x, t, tol, maxit, pltvars);
1105 %! # out = TSTdaspk (outstruct, x, t, tol, maxit, pltvars);
1107 %! out = tst_theta_method (outstruct, x, t, tol, maxit, .5, pltvars, [0 0]);
1108 %! # out = TSTthetamethod (outstruct, x, t, tol, maxit, .5, pltvars, [0 0]);
1110 %! out = tst_odepkg (outstruct, x, t, tol, maxit, pltvars, "oders", [0, 1]);
1111 %! # out = TSTodepkg (outstruct, x, t, tol, maxit, pltvars, "oders", [0, 1]);
1113 %! out = tst_odepkg (outstruct, x, t, tol, maxit, pltvars, "odesx", [0, 1]);
1114 %! # out = TSTodepkg (outstruct, x, t, tol, maxit, pltvars, "odesx", [0, 1]);
1116 %! out = tst_odepkg (outstruct, x, t, tol, maxit, pltvars, "ode2r", [0, 1]);
1117 %! # out = TSTodepkg (outstruct, x, t, tol, maxit, pltvars, "ode2r", [0, 1]);
1119 %! out = tst_odepkg (outstruct, x, t, tol, maxit, pltvars, "ode5r", [0, 1])
1120 %! # out = TSTodepkg (outstruct, x, t, tol, maxit, pltvars, "ode5r", [0, 1])
1122 %! error ("unknown solver");
1125 %! #utl_plot_by_name (t, out, outstruct, pltvars);
1126 %! utl_plot_by_name (t, out, outstruct, pltvars);