home *** CD-ROM | disk | FTP | other *** search
/ NetNews Usenet Archive 1993 #3 / NN_1993_3.iso / spool / alt / sources / 3103 < prev    next >
Encoding:
Text File  |  1993-01-28  |  61.2 KB  |  2,113 lines
  1. Path: sparky!uunet!crdgw1!rpi!usc!howland.reston.ans.net!bogus.sura.net!darwin.sura.net!ukma!cs.widener.edu!dsinc!ub!galileo.cc.rochester.edu!ee.rochester.edu!rbc!al
  2. From: al@rbc.uucp (Al Davis)
  3. Newsgroups: alt.sources
  4. Subject: ACS  circuit simulator  part 18/20
  5. Message-ID: <1993Jan27.041057.12044@rbc.uucp>
  6. Date: 27 Jan 93 04:10:57 GMT
  7. Sender: al@rbc.uucp (Al Davis)
  8. Organization: Huh?
  9. Lines: 2102
  10.  
  11. #! /bin/sh
  12. # This is a shell archive, meaning:
  13. # 1. Remove everything above the #! /bin/sh line.
  14. # 2. Save the resulting text in a file.
  15. # 3. Execute the file with /bin/sh (not csh) to create the files:
  16. #    src/tr_clear.c
  17. #    src/tr_conck.c
  18. #    src/tr_fill.c
  19. #    src/tr_fix.c
  20. #    src/tr_load.c
  21. #    src/tr_out.c
  22. #    src/tr_probe.c
  23. #    src/tr_reviw.c
  24. #    src/tr_setup.c
  25. #    src/tr_solve.c
  26. #    src/tr_sweep.c
  27. # This archive created: Tue Jan 26 22:51:12 1993
  28. export PATH; PATH=/bin:$PATH
  29. if test -f 'src/tr_clear.c'
  30. then
  31.     echo shar: will not over-write existing file "'src/tr_clear.c'"
  32. else
  33. cat << \SHAR_EOF > 'src/tr_clear.c'
  34. /* tr_clear  12/30/92
  35.  * Copyright 1983-1992   Albert Davis
  36.  * Clears working arrays allocated by tralloc.
  37.  * Saves the old time and last iteration.  It is necessary to clear them
  38.  * separately because they must be cleared for each new input, and the old
  39.  * constants saved, and not reallocated.  We cannot use free/calloc because it
  40.  * would require setting new pointer tables.
  41.  *
  42.  * Also, a bunch of other stuff that is here only because of shared variables.
  43.  */
  44. #include "ecah.h"
  45. #include "branch.h"
  46. #include "error.h"
  47. #include "mode.h"
  48. #include "status.h"
  49. #include "declare.h"
  50. /*--------------------------------------------------------------------------*/
  51.     void tr_clear(int);
  52.     void trkeep(void);
  53.     void trestor(void);
  54.     void cmd_keep(void);
  55.     void cmd_unkeep(void);
  56. /*--------------------------------------------------------------------------*/
  57. extern double last_time;        /* time corresp to last voltage        */
  58. extern double *volts;            /* last voltage vector            */
  59. extern double *recon,*imcon;
  60. extern double *oldcon;            /* constant terms            */
  61. extern double *reals, *imags;        /* big array allocation base        */
  62. extern int decomp;            /* how far y-matrix decomposed        */
  63. extern const struct status stats;
  64. extern const unsigned aspace;        /* non-zero array elements        */
  65. extern const double trtime;        /* time at this iteration        */
  66. extern const int inc_mode;        /* make incremental changes        */
  67. extern const int sim_mode;        /* what type simulation is this        */
  68. static int keepold;            /* flag: keep old data            */
  69. /*--------------------------------------------------------------------------*/
  70. void tr_clear(amt)
  71. int amt;
  72. {
  73.  static double last_iter_time;
  74.  
  75.  if (!inc_mode){            /* Clear working array            */
  76.     (void)memset((void*)reals, 0, (size_t)aspace * sizeof(double));
  77.     decomp = amt;
  78.  }
  79.  if (sim_mode == sTRAN  &&  stats.iter[iSTEP] == 1){
  80.     if (trtime > last_iter_time){
  81.        int ii;
  82.        for (ii = 1;  ii <= stats.total_nodes;  ii++){
  83.           oldcon[ii] = imcon[ii];
  84.        }
  85.        last_iter_time = trtime;
  86.     }else{
  87.        /* don't save voltages.  They're wrong !*/
  88.        last_iter_time = trtime;
  89.     }
  90.  }
  91.  if (!inc_mode){            /* Clear new right side            */
  92.     int ii;
  93.     for (ii = 1;  ii <= stats.total_nodes;  ii++){
  94.        recon[ii] = 0.;
  95.     }
  96.  }
  97. }
  98. /*--------------------------------------------------------------------------*/
  99. void trkeep()
  100. {
  101.  if (!keepold){
  102.     int ii;
  103.     for (ii = 1;  ii <= stats.total_nodes;  ii++){
  104.        volts[ii] = imcon[ii];
  105.     }
  106.     last_time = (trtime>0.) ? trtime : 0.;
  107.  }
  108. }
  109. /*--------------------------------------------------------------------------*/
  110. void trestor()
  111. {
  112.  int ii;
  113.  for (ii = 1;  ii <= stats.total_nodes;  ii++){
  114.     imcon[ii] = volts[ii];
  115.  }
  116. }
  117. /*--------------------------------------------------------------------------*/
  118. void cmd_keep()
  119. {
  120.  keepold = YES;
  121. }
  122. /*--------------------------------------------------------------------------*/
  123. void cmd_unkeep()
  124. {
  125.  keepold= NO;
  126. }
  127. /*--------------------------------------------------------------------------*/
  128. /*--------------------------------------------------------------------------*/
  129. SHAR_EOF
  130. fi # end of overwriting check
  131. if test -f 'src/tr_conck.c'
  132. then
  133.     echo shar: will not over-write existing file "'src/tr_conck.c'"
  134. else
  135. cat << \SHAR_EOF > 'src/tr_conck.c'
  136. /* tr_conck.c  12/08/91
  137.  * Copyright 1983-1992   Albert Davis
  138.  * Check convergence, and other stuff done after evaluation
  139.  * to evaluate performance
  140.  */
  141. #include "ecah.h"
  142. #include "branch.h"
  143. #include "convstat.h"
  144. #include "error.h"
  145. #include "io.h"
  146. #include "options.h"
  147. #include "declare.h"
  148. /*--------------------------------------------------------------------------*/
  149.     int    conv_check(const branch_t*);
  150.     int     conchk(double,double,double);
  151. /*--------------------------------------------------------------------------*/
  152. extern const struct ioctrl io;
  153. extern const struct options opt;
  154. /*--------------------------------------------------------------------------*/
  155. /* conv_check: check branch for convergence
  156.  * should be a macro, for speed
  157.  */
  158. int conv_check(brh)
  159. const branch_t *brh;
  160. {
  161.  int conv;
  162.    
  163.     conv = conchk(brh->y1.f0,brh->y0.f0,opt.abstol);
  164.     if ((!io.suppresserrors) && (conv & cNONCONVERGE)){
  165.        error(bPICKY,"%s: non-convergence: (f0) (%s,%s,%s)\n",
  166.           printlabel(brh,NO),
  167.           ftos(brh->y0.f0, "", 7, 0),
  168.           ftos(brh->y1.f0, "", 7, 0),
  169.           ftos(brh->y2.f0, "", 7, 0));
  170.        if (conchk(brh->y1.f1,brh->y0.f1,opt.abstol) & cNONCONVERGE){
  171.           error(bPICKY,"%s: non-convergence: (f1) (%s,%s,%s)\n",
  172.          printlabel(brh,NO),
  173.          ftos(brh->y0.f1, "", 7, 0),
  174.              ftos(brh->y1.f1, "", 7, 0),
  175.          ftos(brh->y2.f1, "", 7, 0));
  176.        }else{
  177.           error(bPICKY,"%s: f1 is ok\n", printlabel(brh,NO));
  178.        }
  179.     }
  180.  
  181.  
  182.  if (!(conv & cNONCONVERGE)){
  183.     conv = conchk(brh->y1.f1,brh->y0.f1,opt.abstol);
  184.     if ((!io.suppresserrors) && (conv & cNONCONVERGE)){
  185.        error(bPICKY,"%s: non-convergence: (f1 only) (%s,%s,%s)\n",
  186.       printlabel(brh,NO),
  187.       ftos(brh->y0.f1, "", 7, 0),
  188.           ftos(brh->y1.f1, "", 7, 0),
  189.       ftos(brh->y2.f1, "", 7, 0));
  190.     }
  191.  }
  192.  return conv;
  193. }
  194. /*--------------------------------------------------------------------------*/
  195. /* conchk: check for convergence
  196.  * should be macro.
  197.  */
  198. int conchk(old,new,abstol)
  199. double old,new,abstol;
  200. {
  201.  if (FABS(new-old) <= abstol)
  202.     return cGOOD;
  203.  else if (FABS(old) > abstol  &&  inorder(opt.lowlim, new/old, opt.uplim))
  204.     return cGOOD;
  205.  else
  206.     return cNONCONVERGE;
  207. }
  208. /*--------------------------------------------------------------------------*/
  209. /*--------------------------------------------------------------------------*/
  210. SHAR_EOF
  211. fi # end of overwriting check
  212. if test -f 'src/tr_fill.c'
  213. then
  214.     echo shar: will not over-write existing file "'src/tr_fill.c'"
  215. else
  216. cat << \SHAR_EOF > 'src/tr_fill.c'
  217. /* tr_fill  03/26/92
  218.  * Copyright 1983-1992   Albert Davis
  219.  * Fill working matrix for transient and dc analysis
  220.  * Returns convergence status.
  221.  */
  222. #include "ecah.h"
  223. #include "branch.h"
  224. #include "convstat.h"
  225. #include "dev.h"
  226. #include "mode.h"
  227. #include "status.h"
  228. #include "declare.h"
  229. /*--------------------------------------------------------------------------*/
  230.     int    tr_fill(void);
  231.     int    tr_fill_rl(branch_t*);
  232.     int    tr_fill_ll(branch_t**);
  233.     void    tr_unfill_rl(branch_t*);
  234. /*--------------------------------------------------------------------------*/
  235. extern struct status stats;
  236. /*--------------------------------------------------------------------------*/
  237. /* tr_fill:  process netlist, evaluate models, fill matrix, etc
  238.  * clears matrix, re-fills it, checks for convergence
  239.  * does not evaluate matrix
  240.  */
  241. int tr_fill()
  242. {
  243.  int conv;                /* convergence flag                */
  244.  
  245.  time_start(&(stats.load));
  246.  conv = tr_fill_rl(firstbranch_dev());
  247.  if (stats.iter[iSTEP] == 1)
  248.     conv = cNONCONVERGE;
  249.  time_stop(&(stats.load));
  250.  return conv;
  251. }
  252. /*--------------------------------------------------------------------------*/
  253. /* tr_fill_rl: evaluate a list of models, fill matrix
  254.  * evaluates a list (or sublist), checks convergence, etc.
  255.  * argument is the head of the netlist.
  256.  * recursively called to evaluate subcircuits
  257.  */
  258. int tr_fill_rl(brh)
  259. branch_t *brh;
  260. {
  261.  int conv;
  262.  const branch_t *stop;
  263.  
  264.  conv = cGOOD;
  265.  if (exists(brh)){
  266.     stop = brh;
  267.     do {
  268.        conv |= dotr_branch(brh);
  269.     } while (brh=nextbranch_dev(brh),  brh != stop);
  270.  }
  271.  return conv;
  272. }
  273. /*--------------------------------------------------------------------------*/
  274. /* tr_fill_ll: linear list version of tr_fill_rl
  275.  */
  276. int tr_fill_ll(bl)
  277. branch_t **bl;
  278. {
  279.  int conv;
  280.  branch_t *brh;
  281.  
  282.  conv = cGOOD;
  283.  if (bl){
  284.     for (brh = *bl;  exists(brh);  ++bl){
  285.        conv |= dotr_branch(brh);
  286.     }
  287.  }
  288.  return conv;
  289. }
  290. /*--------------------------------------------------------------------------*/
  291. /* tr_unfill_rl: remove evaluate a list of models from the matrix
  292.  * recursively called to unevaluate subcircuits
  293.  */
  294. void tr_unfill_rl(brh)
  295. branch_t *brh;
  296. {
  297.  const branch_t *stop;
  298.  if (exists(brh)){
  299.     stop = brh;
  300.     do {
  301.        untr_branch(brh);
  302.     } while (brh=nextbranch_dev(brh),  brh != stop);
  303.  }
  304. }
  305. /*--------------------------------------------------------------------------*/
  306. /*--------------------------------------------------------------------------*/
  307. SHAR_EOF
  308. fi # end of overwriting check
  309. if test -f 'src/tr_fix.c'
  310. then
  311.     echo shar: will not over-write existing file "'src/tr_fix.c'"
  312. else
  313. cat << \SHAR_EOF > 'src/tr_fix.c'
  314. /* trfix  01/18/93
  315.  * Copyright 1983-1992   Albert Davis
  316.  * Takes care of nonlinearities, behavioral modeling, etc
  317.  *     in transient and dc analysis.
  318.  */
  319. #include "ecah.h"
  320. #include "branch.h"
  321. #include "error.h"
  322. #include "expr.h"
  323. #include "declare.h"
  324. /*--------------------------------------------------------------------------*/
  325.     void    trfix0(branch_t*);
  326.     void    trfix1(branch_t*);
  327. static    void    trf_ac(const branch_t*,double**,int*,cpoly1_t*);
  328. static    void    trf_dc(const branch_t*,double**,int*,cpoly1_t*);
  329. static    void    trf_dctran(const branch_t*,double**,int*,cpoly1_t*);
  330. static    void    trf_frequency(const branch_t*,double**,int*,cpoly1_t*);
  331. static    void    trf_period(const branch_t*,double**,int*,cpoly1_t*);
  332. static    void    trf_ramp(const branch_t*,double**,int*,cpoly1_t*);    
  333. static    void    trf_time(const branch_t*,double**,int*,cpoly1_t*);    
  334. static    void    trf_tran(const branch_t*,double**,int*,cpoly1_t*);    
  335. static    void    trf_bandwidth(const branch_t*,double**,cpoly1_t*);
  336. static    void    trf_complex(const branch_t*,double**,cpoly1_t*);
  337. static    void    trf_cornerdown(const branch_t*,double**,cpoly1_t*);
  338. static    void    trf_cornerup(const branch_t*,double**,cpoly1_t*);
  339. static    void    trf_delay(const branch_t*,double**,cpoly1_t*);    
  340. static    void    trf_exp(const branch_t*,double**,cpoly1_t*);    
  341. static    void    trf_expterm(const branch_t*,double**,cpoly1_t*);
  342. static    void    trf_generator(const branch_t*,double**,cpoly1_t*);
  343. static    void    trf_max(const branch_t*,double**,cpoly1_t*);    
  344. static    void    trf_netfunc(const branch_t*,double**,cpoly1_t*);    
  345. static    void    trf_notch(const branch_t*,double**,cpoly1_t*);    
  346. static    void    trf_numeric(const branch_t*,double**,cpoly1_t*);    
  347. static    void    trf_offset(const branch_t*,double**,cpoly1_t*);
  348. static    void    trf_polar(const branch_t*,double**,cpoly1_t*);    
  349. static    void    trf_polyterm(const branch_t*,double**,cpoly1_t*);
  350. static    void    trf_pulse(const branch_t*,double**,cpoly1_t*);    
  351. static    void    trf_pwl(const branch_t*,double**,cpoly1_t*);    
  352. static    void    trf_sffm(const branch_t*,double**,cpoly1_t*);
  353. static    void    trf_sin(const branch_t*,double**,cpoly1_t*);    
  354. static    void    trf_tanh(const branch_t*,double**,cpoly1_t*);    
  355. static    void    trf_ic(const branch_t*,double**,cpoly1_t*);    
  356. static    void    trf_ii(const branch_t*,double**,cpoly1_t*);    
  357. static    void    trf_iv(const branch_t*,double**,cpoly1_t*);    
  358. static    void    trf_tempco(const branch_t*,double**,cpoly1_t*);
  359. /*--------------------------------------------------------------------------*/
  360. #define arg0 (arg[0][0])
  361. #define arg1 (arg[0][1])
  362. #define arg2 (arg[0][2])
  363. #define arg3 (arg[0][3])
  364. #define arg4 (arg[0][4])
  365. #define arg5 (arg[0][5])
  366. #define arg6 (arg[0][6])
  367. #define arg7 (arg[0][7])
  368.  
  369. extern const double genout;    /* the bench signal generator            */
  370. extern const double trtime;    /* transient analysis time            */
  371.  
  372. extern int in_curr;        /* flag: initial current            */
  373. extern int in_volt;        /* flag: initial voltage            */
  374. extern int in_cond;        /* flag: generic initial condition        */
  375. extern double init_curr;    /* initial current                */
  376. extern double init_volt;    /* initial voltage                */
  377. extern double init_cond;    /* generic initial condition            */
  378.  
  379. static int skip;
  380. static double periodtime;    /* time since beginning of period        */
  381. static double reltime;        /* time since event                */
  382. /*--------------------------------------------------------------------------*/
  383. void trfix0(brh)
  384. branch_t *brh;
  385. {
  386.  brh->y0.f0 = NOT_VALID;
  387. }
  388. /*--------------------------------------------------------------------------*/
  389. void trfix1(brh)
  390. branch_t *brh;
  391. {
  392.  cpoly1_t y;            /* return value                    */
  393.                 /* y.c1 is the slope, y.c0 is y intercept   */
  394.  skip = NO;
  395.  periodtime = reltime = trtime;
  396.  y.x  = brh->y0.x;
  397.  y.f1 = brh->val;
  398.  y.c0 = 0.;
  399.  
  400.  if (brh->x){
  401.     struct expr *x;
  402.     double *arg;
  403.     int *key;
  404.  
  405.     x = (struct expr*)brh->x;
  406.     arg = x->args->args;
  407.     for (key = x->keys->args;  *key;  key++){
  408.        switch (*key){
  409.       case eAC:        trf_ac(brh,&arg,&skip,&y);     break;
  410.       case eDC:        trf_dc(brh,&arg,&skip,&y);     break;
  411.       case eDCTRAN:        trf_dctran(brh,&arg,&skip,&y);     break;
  412.       case eFREQUENCY:    trf_frequency(brh,&arg,&skip,&y);break;
  413.       case ePERIOD:        trf_period(brh,&arg,&skip,&y);     break;
  414.       case eRAMP:        trf_ramp(brh,&arg,&skip,&y);     break;
  415.       case eTIME:        trf_time(brh,&arg,&skip,&y);     break;
  416.       case eTRAN:        trf_tran(brh,&arg,&skip,&y);     break;
  417.  
  418.       case eBANDWIDTH:    trf_bandwidth(brh,&arg,&y);     break;
  419.       case eCOMPLEX:    trf_complex(brh,&arg,&y);     break;
  420.       case eCORNERDOWN:    trf_cornerdown(brh,&arg,&y);     break;
  421.       case eCORNERUP:    trf_cornerup(brh,&arg,&y);     break;
  422.       case eDELAY:        trf_delay(brh,&arg,&y);         break;
  423.       case eEXP:        trf_exp(brh,&arg,&y);         break;
  424.       case eEXPTERM:    trf_expterm(brh,&arg,&y);     break;
  425.       case eGENERATOR:    trf_generator(brh,&arg,&y);     break;
  426.       case eMAX:        trf_max(brh,&arg,&y);         break;
  427.       case eNETFUNC:    trf_netfunc(brh,&arg,&y);     break;
  428.       case eNOTCH:        trf_notch(brh,&arg,&y);         break;
  429.       case eNUMERIC:    trf_numeric(brh,&arg,&y);     break;
  430.       case eOFFSET:        trf_offset(brh,&arg,&y);     break;
  431.       case ePOLAR:        trf_polar(brh,&arg,&y);         break;
  432.       case ePOLYTERM:    trf_polyterm(brh,&arg,&y);     break;
  433.       case ePULSE:        trf_pulse(brh,&arg,&y);         break;
  434.       case ePWL:        trf_pwl(brh,&arg,&y);         break;
  435.       case eSFFM:        trf_sffm(brh,&arg,&y);         break;
  436.       case eSIN:        trf_sin(brh,&arg,&y);         break;
  437.       case eTANH:        trf_tanh(brh,&arg,&y);         break;
  438.  
  439.       case eIC:        trf_ic(brh,&arg,&y);         break;
  440.       case eII:        trf_ii(brh,&arg,&y);         break;
  441.       case eIV:        trf_iv(brh,&arg,&y);         break;
  442.       case eTEMPCO:        trf_tempco(brh,&arg,&y);     break;
  443.  
  444.       default:  error(bWARNING, "%s: undefined function: %d\n", printlabel(brh,NO), *key); break;
  445.        }
  446.     }
  447.  }
  448.  brh->y0.x  = y.x;
  449.  brh->y0.f0 = y.c0 + y.x * y.f1;
  450.  brh->y0.f1 = y.f1;
  451. }
  452. /*--------------------------------------------------------------------------*/
  453. /* trf_ac:  keyword: ac
  454.  * following args for ac analysis only
  455.  * here : skip them
  456.  */
  457. /*ARGSUSED*/
  458. static void trf_ac(brh,arg,skip,y)        /* no args        */
  459. const branch_t *brh;
  460. double **arg;
  461. int *skip;
  462. cpoly1_t* y;
  463. {
  464.  *skip = YES;
  465.  *arg += aAC;
  466. }
  467. /*--------------------------------------------------------------------------*/
  468. /* trf_dc:  keyword: dc
  469.  * following args for dc and transient analysis
  470.  */
  471. /*ARGSUSED*/
  472. static void trf_dc(brh,arg,skip,y)        /* no args        */
  473. const branch_t *brh;
  474. double **arg;
  475. int *skip;
  476. cpoly1_t *y;
  477. {
  478.  y->f1 = y->c0 = 0.;
  479.  *skip = NO;
  480.  *arg += aDC;
  481. }
  482. /*--------------------------------------------------------------------------*/
  483. /* trf_dctran:  keyword: dctran
  484.  * following args for dc and transient analysis
  485.  * this does apply to initial conditions.
  486.  */
  487. /*ARGSUSED*/
  488. static void trf_dctran(brh,arg,skip,y)        /* no args        */
  489. const branch_t *brh;
  490. double **arg;
  491. int *skip;
  492. cpoly1_t* y;
  493. {
  494.  y->f1 = y->c0 = 0.;
  495.  *skip = NO;
  496.  *arg += aDCTRAN;
  497. }
  498. /*--------------------------------------------------------------------------*/
  499. /* trf_frequency:  keyword: frequency
  500.  * the value is frequency dependent (ac only)
  501.  * works only in ac analysis, otherwise could be non-causal.
  502.  */
  503. /*ARGSUSED*/
  504. static void trf_frequency(brh,arg,skip,y)    /* arg0 = test freq    */
  505. const branch_t *brh;
  506. double **arg;
  507. int *skip;
  508. cpoly1_t *y;
  509. {
  510.  *skip = YES;
  511.  *arg += aFREQUENCY;
  512. }
  513. /*--------------------------------------------------------------------------*/
  514. /* trf_period:  keyword: period
  515.  * the component is periodic in time
  516.  * periodtime and reltime are reset every period
  517.  * all time dependencies are based on reltime and periodtime
  518.  * for sensible behavior, this key should be first.
  519.  * implies that following args for transient analysis only.
  520.  */
  521. /*ARGSUSED*/
  522. static void trf_period(brh,arg,skip,y)        /* arg0 = period    */
  523. const branch_t *brh;
  524. double **arg;
  525. int *skip;
  526. cpoly1_t *y;
  527. {
  528.  if (trtime >= 0.){ /* if transient analysis */
  529.     y->f1 = y->c0 = 0.;
  530.     *skip = NO;
  531.     periodtime = reltime = (arg0 == 0.)  ?  trtime  :  fmod(trtime, arg0);
  532.  }else{
  533.     *skip = YES;
  534.  }
  535.  *arg += aPERIOD;
  536. }
  537. /*--------------------------------------------------------------------------*/
  538. /*ARGSUSED*/
  539. static void trf_ramp(brh,arg,skip,y)
  540. const branch_t *brh;
  541. double **arg;
  542. int *skip;
  543. cpoly1_t *y;
  544. {
  545.  error(bWARNING,"ramp not implemented: %s\n", printlabel(brh,NO));
  546.  *arg += aRAMP;
  547. }
  548. /*--------------------------------------------------------------------------*/
  549. /* trf_time:  keyword: time
  550.  * following args take effect at the given time (switch)
  551.  * obviously, transient only.
  552.  */
  553. /*ARGSUSED*/
  554. static void trf_time(brh,arg,skip,y)        /* arg0 = switch time    */
  555. const branch_t *brh;
  556. double **arg;
  557. int *skip;
  558. cpoly1_t *y;
  559. {
  560.  if (periodtime >= arg0){
  561.     y->f1 = y->c0 = 0.;
  562.     reltime = periodtime - arg0;
  563.     *skip = NO;
  564.  }else{
  565.     *skip = YES;
  566.  }
  567.  *arg += aTIME;
  568. }
  569. /*--------------------------------------------------------------------------*/
  570. /* trf_tran:  keyword: transient
  571.  * following args for transient analysis only.
  572.  * applies to the dc analysis that computes initial conditions,
  573.  * but not the steady state dc analysis.
  574.  */
  575. /*ARGSUSED*/
  576. static void trf_tran(brh,arg,skip,y)        /* no args        */
  577. const branch_t *brh;
  578. double **arg;
  579. int *skip;
  580. cpoly1_t *y;
  581. {
  582.  if (trtime >= 0.){ /* if transient analysis */
  583.     y->f1 = y->c0 = 0.;
  584.     *skip = NO;
  585.  }else{
  586.     *skip = YES;
  587.  }
  588.  *arg += aTRAN;
  589. }
  590. /*--------------------------------------------------------------------------*/
  591. /* trf_bandwidth:  function: bandwidth
  592.  * gain block, with a bandwidth. (ac only)
  593.  * bad design: should be keyword instead.
  594.  * not implemented in transient analysis
  595.  * but could be as R-C or similar
  596.  */
  597. static void trf_bandwidth(brh,arg,y)        /* arg0 = dc gain    */
  598. const branch_t *brh;                /* arg1 = 3 db freq.    */
  599. double **arg;
  600. cpoly1_t *y;
  601. {
  602.  if (!skip){
  603.     y->f1 += arg0;
  604.     if (trtime > 0.)
  605.     error(bPICKY,"%s: bandwidth not supported in transient analysis\n",printlabel(brh,NO));
  606.  }
  607.  *arg += aBANDWIDTH;
  608. }
  609. /*--------------------------------------------------------------------------*/
  610. /* trf_complex:  function: complex
  611.  * complex value, cartesian coordinates, in frequency domain (ac only)
  612.  * not possible in transient or dc because would be non-causal
  613.  */
  614. static void trf_complex(brh,arg,y)        /* arg0 = real part     */
  615. const branch_t *brh;                /* arg1 = imaginary part */
  616. double **arg;
  617. cpoly1_t *y;
  618. {
  619.  if (!skip){
  620.     y->f1 += arg0;
  621.     error(bPICKY, "non-causal circuit: %s\n", printlabel(brh,NO));
  622.  }
  623.  *arg += aCOMPLEX;
  624. }
  625. /*--------------------------------------------------------------------------*/
  626. /* trf_cornerdown:  function: cornerdown
  627.  * piecewise linear function:
  628.  * output is 0 for in >= arg1
  629.  * derivative is arg0 for in <= arg1
  630.  */
  631. static void trf_cornerdown(brh,arg,y)        /* arg0 = active slope    */
  632. const branch_t *brh;                /* arg1 = break point    */
  633. double **arg;
  634. cpoly1_t *y;
  635. {
  636.  if (!skip  &&  y->x < arg1){
  637.     y->f1 += arg0;
  638.     y->c0 -= arg1 * arg0;
  639.  }
  640.  *arg += aCORNERDOWN;
  641. }
  642. /*--------------------------------------------------------------------------*/
  643. /* trf_cornerup:  function: cornerup
  644.  * piecewise linear function:
  645.  * output is 0 for in <= arg1
  646.  * derivative is arg0 for in >= arg1
  647.  */
  648. static void trf_cornerup(brh,arg,y)        /* arg0 = active slope    */
  649. const branch_t *brh;                /* arg1 = break point    */
  650. double **arg;
  651. cpoly1_t *y;
  652. {
  653.  if (!skip  &&  y->x > arg1){
  654.     y->f1 += arg0;
  655.     y->c0 -= arg1 * arg0;
  656.  }
  657.  *arg += aCORNERUP;
  658. }
  659. /*--------------------------------------------------------------------------*/
  660. /* trf_delay:  function: delay
  661.  * time delay (ac only)
  662.  * not implemented, needs storage of past time values
  663.  */
  664. static void trf_delay(brh,arg,y)        /* arg0 = magnitude    */
  665. const branch_t *brh;                /* arg1 = delay        */
  666. double **arg;
  667. cpoly1_t *y;
  668. {
  669.  if (!skip){
  670.     y->f1 += arg0;
  671.     error(bPICKY,"%s: delay not supported in transient or analysis\n",printlabel(brh,NO));
  672.  }
  673.  *arg += aDELAY;
  674. }
  675. /*--------------------------------------------------------------------------*/
  676. /* trf_exp:  spice compatible function: exp
  677.  * spice source exponential function: exponential function of time
  678.  * known bugs: defaults wrong: could divide by zero
  679.  */
  680. /*ARGSUSED*/
  681. static void trf_exp(brh,arg,y)            /* arg0 = initial value      */
  682. const branch_t *brh;                /* arg1 = pulsed value      */
  683. double **arg;                    /* arg2 = rise delay      */
  684. cpoly1_t *y;                    /* arg3 = rise time const */
  685. {                        /* arg4 = fall delay      */
  686.  if (!skip){                    /* arg5 = fall time const */
  687.     y->f1 += arg0;
  688.     if (reltime > arg2)
  689.        y->f1 += (arg1 - arg0) * (1. - exp(-(reltime-arg2)/arg3));
  690.     if (reltime > arg4)
  691.        y->f1 += (arg0 - arg1) * (1. - exp(-(reltime-arg4)/arg5));
  692.  }
  693.  *arg += aEXP;
  694. }
  695. /*--------------------------------------------------------------------------*/
  696. /* trf_expterm:  function: expterm
  697.  * exponent term, non-integer power term (not exponential)
  698.  * like polyterm, but for non-integers
  699.  * only works for positive input, because negative number raised to
  700.  * non-integer power is usually complex.
  701.  */
  702. static void trf_expterm(brh,arg,y)        /* arg0 = coefficient    */
  703. const branch_t *brh;                /* arg1 = exponent    */
  704. double **arg;
  705. cpoly1_t *y;
  706. {
  707.  if (!skip  &&  y->x > 0.){
  708.     double coeff;
  709.     coeff = arg0 * pow(y->x,arg1-1);
  710.     y->f1 += coeff * arg1;
  711.     y->c0 += coeff * y->x * (1-arg1);
  712.  }
  713.  *arg += aEXPTERM;
  714. }
  715. /*--------------------------------------------------------------------------*/
  716. /* trf_generator:  function: generator
  717.  * value is derived from the "signal generator" (generator command)
  718.  * intended for fixed sources, as circuit input.
  719.  * component specific period determines only whether to use this or not
  720.  * it does not actually affect the operation of the signal generator
  721.  */
  722. /*ARGSUSED*/
  723. static void trf_generator(brh,arg,y)        /* arg0 = scale factor    */
  724. const branch_t *brh;
  725. double **arg;
  726. cpoly1_t *y;
  727. {
  728.  if (!skip)
  729.     y->f1 += arg0 * genout;
  730.  *arg += aGENERATOR;
  731. }
  732. /*--------------------------------------------------------------------------*/
  733. /* trf_max:  function: max
  734.  * piecewise linear function: block that clips
  735.  * bad design: should be keyword, so it can apply to the total part
  736.  *           should be able to specify upper and lower limits
  737.  */
  738. static void trf_max(brh,arg,y)            /* arg0 = normal value     */
  739. const branch_t *brh;                /* arg1 = output clip pt */
  740. double **arg;
  741. cpoly1_t *y;
  742. {
  743.  if (!skip  &&  arg0 != 0.){
  744.     double clip_input;
  745.     clip_input = arg1/arg0;
  746.     if (fabs(y->x-brh->y1.x) < fabs(y->x-brh->y2.x)*10.){    /* normal */
  747.        if (y->x < -clip_input)                    /* -clip */
  748.       y->c0 -= arg1;
  749.        else if (y->x > clip_input)                /* +clip */
  750.       y->c0 += arg1;
  751.        else                            /* linear */
  752.       y->f1 += arg0;
  753.     }else{                    /* oscillating */
  754.                             /* don't skip middle region */
  755.        if (y->x < -clip_input  &&  brh->y1.x < clip_input)
  756.       y->c0 -= arg1;
  757.        else if (y->x > clip_input  &&  brh->y1.x > -clip_input)
  758.       y->c0 += arg1;
  759.        else
  760.       y->f1 += arg0;
  761.     }
  762.  }
  763.  *arg += aMAX;
  764. }
  765. /*--------------------------------------------------------------------------*/
  766. /* trf_netfunc:  function: netfunction
  767.  * gain block, nth order response. (ac only)
  768.  * not implemented in transient analysis
  769.  * but could be as R-L-C or similar
  770.  */                        /* arg0 = arg count    */
  771. static void trf_netfunc(brh,arg,y)        /* arg1 = dc gain    */
  772. const branch_t *brh;                /* arg3 = order of denom*/
  773. double **arg;                    /* arg* = coefs of denom*/
  774. cpoly1_t *y;                    /* arg* = coefs of numer*/
  775. {
  776.  int argcount;
  777.  argcount = (int)arg0;
  778.  
  779.  if (!skip){
  780.     y->f1 += arg1;
  781.     if (trtime > 0.)
  782.         error(bPICKY,"%s: netfunction not supported in transient analysis\n",printlabel(brh,NO));
  783.  }
  784.  *arg += argcount;
  785. }
  786. /*--------------------------------------------------------------------------*/
  787. /* trf_notch:  function: notch
  788.  * piecewise linear function: crossover notch, dead zone
  789.  * symmetric around zero
  790.  */
  791. static void trf_notch(brh,arg,y)        /* arg0 = normal value     */
  792. const branch_t *brh;                /* arg1 = dead zone size */
  793. double **arg;
  794. cpoly1_t *y;
  795. {
  796.  if (!skip){
  797.     if (y->x > arg1){
  798.        y->f1 += arg0;
  799.        y->c0 -= arg1 * arg0;
  800.     }else if (y->x < -arg1){
  801.        y->f1 += arg0;
  802.        y->c0 += arg1 * arg0;
  803.     }
  804.  }
  805.  *arg += aNOTCH;
  806. }
  807. /*--------------------------------------------------------------------------*/
  808. /* trf_numeric:  simple numeric argument
  809.  */
  810. /*ARGSUSED*/
  811. static void trf_numeric(brh,arg,y)        /* arg0 = value        */
  812. const branch_t *brh;
  813. double **arg;
  814. cpoly1_t *y;
  815. {
  816.  if (!skip)
  817.     y->f1 += arg0;
  818.  *arg += aNUMERIC;
  819. }
  820. /*--------------------------------------------------------------------------*/
  821. /* trf_offset:  function: offset
  822.  * fixed dc offset
  823.  * bad design: should be keyword
  824.  */
  825. /*ARGSUSED*/
  826. static void trf_offset(brh,arg,y)        /* arg0 = gain        */
  827. const branch_t *brh;                /* arg1 = output offset    */
  828. double **arg;
  829. cpoly1_t *y;
  830. {
  831.  if (!skip){
  832.     y->f1 += arg0;
  833.     y->c0 += arg1 * arg0;
  834.  }
  835.  *arg += aOFFSET;
  836. }
  837. /*--------------------------------------------------------------------------*/
  838. /* trf_polar:  function: polar
  839.  * complex value in polar coordinates, in frequency domain (ac only)
  840.  * not possible in transient or dc because would be non-causal
  841.  */
  842. static void trf_polar(brh,arg,y)        /* arg0 = magnitude    */
  843. const branch_t *brh;                /* arg1 = phase        */
  844. double **arg;
  845. cpoly1_t *y;
  846. {
  847.  if (!skip){
  848.     y->f1 += arg0;
  849.     error(bPICKY, "non-causal circuit: %s\n", printlabel(brh,NO));
  850.  }
  851.  *arg += aPOLAR;
  852. }
  853. /*--------------------------------------------------------------------------*/
  854. /* trf_polyterm:  function: polyterm
  855.  * polynomial term, one term in a polynomial
  856.  * power must be integer, is rounded to nearest integer
  857.  * caution: will divide by zero with zero input and negative exponent
  858.  */
  859. static void trf_polyterm(brh,arg,y)             /* arg0 = coefficient   */
  860. const branch_t *brh;                              /* arg1 = exponent      */
  861. double **arg;
  862. cpoly1_t *y;
  863. {
  864.  if (!skip){
  865.     int expo;
  866.     expo = (int)floor(arg1+.5);
  867.     if (expo == 0){
  868.        y->c0 += arg0;
  869.     }else{ /* (expo != 0) */
  870.        double coeff;
  871.        coeff = arg0 * ipow(y->x,expo-1);
  872.        y->f1 += coeff * expo;
  873.        y->c0 += coeff * y->x * (1-expo);
  874.     }
  875.  }
  876.  *arg += aPOLYTERM;
  877. }
  878. /*--------------------------------------------------------------------------*/
  879. /* trf_pulse: spice compatible function: pulse
  880.  * spice pulse function (for sources)
  881.  */
  882. /*ARGSUSED*/
  883. static void trf_pulse(brh,arg,y)        /* arg0 = initial value    */
  884. const branch_t *brh;                /* arg1 = pulsed value    */
  885. double **arg;                    /* arg2 = delay time    */
  886. cpoly1_t *y;                    /* arg3 = rise time    */
  887. {                        /* arg4 = fall time    */
  888.                         /* arg5 = pulse width    */
  889.                         /* arg6 = period    */
  890.  if (!skip){
  891.     double pw;
  892.     double time;
  893.     
  894.     pw   = (arg5 == 0.) ? reltime : arg5;
  895.     time = reltime;
  896.     if (arg6 != 0.){
  897.        while (time >= arg2+arg6){
  898.           time -= arg6;
  899.        }
  900.     }
  901.  
  902.     if (time >= arg2+arg3+pw+arg4){        /* past pulse    */
  903.        y->f1 += arg0;
  904.     }else if (time >= arg2+arg3+pw){        /* falling     */
  905.        double interp;
  906.        interp = (time - (arg2+arg3+pw)) / arg4;
  907.        y->f1 += arg1 + interp * (arg0 - arg1);
  908.     }else if (time >= arg2+arg3){        /* pulse val     */
  909.        y->f1 += arg1;
  910.     }else if (time >= arg2){            /* rising     */
  911.        double interp;
  912.        interp = (time - arg2) / arg3;
  913.        y->f1 += arg0 + interp * (arg1 - arg0);
  914.     }else{                    /* init val    */
  915.        y->f1 += arg0;
  916.     }
  917.  }
  918.  *arg += aPULSE;
  919. }
  920. /*--------------------------------------------------------------------------*/
  921. /* trf_pwl:  spice compatible function: pwl
  922.  * "piece-wise linear" point-plotting function of time.
  923.  * BUG: different from SPICE when first point is not zero
  924.  */
  925. /*ARGSUSED*/                    /* arg0 = arg count    */
  926. static void trf_pwl(brh,arg,y)            /* arg1 = time        */
  927. const branch_t *brh;                /* arg2 = value        */
  928. double **arg;                    /* etc.            */
  929. cpoly1_t *y;
  930. {
  931.  int argcount;
  932.  argcount = (int)(arg[0][0]);
  933.  
  934.  if (!skip){
  935.     int i;
  936.     double lowtime;
  937.     double hitime;
  938.     double lowvalue;
  939.     double hivalue;
  940.  
  941.     lowtime = 0.;
  942.     lowvalue = arg[0][2];
  943.     for (i=1;
  944.          i < argcount  &&  reltime > arg[0][i]  &&  lowtime <= arg[0][i];
  945.     /**/ i += 2){
  946.        lowtime = arg[0][i];
  947.        lowvalue = arg[0][i+1];
  948.     }
  949.  
  950.     hitime = arg[0][i];
  951.     hivalue = arg[0][i+1];
  952.     if (i < argcount  &&  hitime > lowtime){
  953.        double ratio;
  954.        ratio = (reltime - lowtime) / (hitime - lowtime);
  955.        y->f1 += lowvalue + (hivalue - lowvalue) * ratio;
  956.     }else{
  957.        y->f1 += lowvalue;
  958.     }
  959.  }
  960.  *arg += argcount;
  961. }
  962. /*--------------------------------------------------------------------------*/
  963. /* trf_sffm:  spice compatible function: sffm
  964.  * single frequency frequency modulation
  965.  */
  966. /*ARGSUSED*/
  967. static void trf_sffm(brh,arg,y)            /* arg0 = dc offset    */
  968. const branch_t *brh;                /* arg1 = amplitude    */
  969. double **arg;                    /* arg2 = carrier freq    */
  970. cpoly1_t *y;                    /* arg3 = mod index    */
  971. {                        /* arg4 = signal freq    */
  972.  if (!skip){
  973.     double mod;
  974.     mod = arg3 * sin(PI2*arg4*reltime);
  975.     y->f1 += arg0 + arg1 * sin(PI2*arg2*reltime + mod);
  976.  }
  977.  *arg += aSFFM;
  978. }
  979. /*--------------------------------------------------------------------------*/
  980. /* trf_sin:  spice compatible function: sin
  981.  * sinusoidal function, mainly for sources
  982.  * value is decaying sinusoidal function of time
  983.  */
  984. /*ARGSUSED*/
  985. static void trf_sin(brh,arg,y)            /* arg0 = offset    */
  986. const branch_t *brh;                /* arg1 = amplitude    */
  987. double **arg;                    /* arg2 = frequency    */
  988. cpoly1_t *y;                    /* arg3 = delay        */
  989. {                        /* arg4 = damping    */
  990.  if (!skip){
  991.     y->f1 += arg0;
  992.     if (reltime > arg3){
  993.        double x;
  994.        x = arg1 * sin(PI2*arg2*(reltime-arg3));
  995.        if (arg4 != 0.)
  996.           x *= exp(-(reltime-arg3)*arg4);
  997.        y->f1 += x;
  998.     }
  999.  }
  1000.  *arg += aSIN;
  1001. }
  1002. /*--------------------------------------------------------------------------*/
  1003. /* trf_tanh:  function: tanh
  1004.  * piecewise linear function: block that clips
  1005.  * bad design: should be keyword, so it can apply to the total part
  1006.  *           should be able to specify upper and lower limits
  1007.  */
  1008. static void trf_tanh(brh,arg,y)            /* arg0 = normal value     */
  1009. const branch_t *brh;                /* arg1 = output clip pt */
  1010. double **arg;
  1011. cpoly1_t *y;
  1012. {
  1013.  if (!skip){
  1014.     if (arg1 != 0.){
  1015.        double aa;
  1016.        double cosine;
  1017.        aa = y->x * arg0/arg1;
  1018.        cosine = cosh(aa);
  1019.        y->f1 += arg0 / (cosine*cosine);
  1020.        y->c0 += arg1 * tanh(aa) - y->f1 * y->x;
  1021.     }
  1022.     /* else 0 */
  1023.  }
  1024.  *arg += aTANH;
  1025. }
  1026. /*--------------------------------------------------------------------------*/
  1027. /*ARGSUSED*/
  1028. static void trf_ic(brh,arg,y)
  1029. const branch_t *brh;
  1030. double **arg;
  1031. cpoly1_t *y;
  1032. {
  1033.  if (!skip){
  1034.     in_cond = YES;
  1035.     init_cond = arg0;
  1036.  }
  1037.  *arg += 1;
  1038. }
  1039. /*--------------------------------------------------------------------------*/
  1040. /*ARGSUSED*/
  1041. static void trf_ii(brh,arg,y)
  1042. const branch_t *brh;
  1043. double **arg;
  1044. cpoly1_t *y;
  1045. {
  1046.  if (!skip){
  1047.     in_curr = YES;
  1048.     init_curr = arg0;
  1049.  }
  1050.  *arg += 1;
  1051. }
  1052. /*--------------------------------------------------------------------------*/
  1053. /*ARGSUSED*/
  1054. static void trf_iv(brh,arg,y)
  1055. const branch_t *brh;
  1056. double **arg;
  1057. cpoly1_t *y;
  1058. {
  1059.  if (!skip){
  1060.     in_volt = YES;
  1061.     init_volt = arg0;
  1062.  }
  1063.  *arg += 1;
  1064. }
  1065. /*--------------------------------------------------------------------------*/
  1066. /*ARGSUSED*/
  1067. static void trf_tempco(brh,arg,y)
  1068. const branch_t *brh;
  1069. double **arg;
  1070. cpoly1_t *y;
  1071. {
  1072. #ifdef NEVER
  1073.  tempco = arg0;
  1074. #endif
  1075.  *arg += 1;
  1076. }
  1077. /*--------------------------------------------------------------------------*/
  1078. /*--------------------------------------------------------------------------*/
  1079. SHAR_EOF
  1080. fi # end of overwriting check
  1081. if test -f 'src/tr_load.c'
  1082. then
  1083.     echo shar: will not over-write existing file "'src/tr_load.c'"
  1084. else
  1085. cat << \SHAR_EOF > 'src/tr_load.c'
  1086. /* tr_load  03/11/92
  1087.  * Copyright 1983-1992   Albert Davis
  1088.  * Load matrix from pre-computed values
  1089.  */
  1090. #include "ecah.h"
  1091. #include "array.h"
  1092. #include "branch.h"
  1093. #include "error.h"
  1094. #include "mode.h"
  1095. #include "options.h"
  1096. #include "status.h"
  1097. #include "declare.h"
  1098. /*--------------------------------------------------------------------------*/
  1099.     void    trloadsource(branch_t*);
  1100.     void    trloadpassive(branch_t*);
  1101.     void    trloadactive(branch_t*);
  1102.     void    unloadsource(branch_t*);
  1103.     void    unloadpassive(branch_t*);
  1104.     void    unloadactive(branch_t*);
  1105. /*--------------------------------------------------------------------------*/
  1106. extern const struct options opt;
  1107. extern       struct status stats;
  1108. extern int inc_mode;            /* make incremental changes            */
  1109. extern double *recon;            /* constant terms                */
  1110. /*--------------------------------------------------------------------------*/
  1111. /* trloadsource: load source to right side vector
  1112.  */
  1113. void trloadsource(brh)
  1114. branch_t *brh;
  1115. {
  1116.  double dc0;
  1117.  if (brh->loaditer == stats.iter[iTOTAL])
  1118.     error(bDANGER, "%s: double load\n", printlabel(brh,0));
  1119.  brh->loaditer = stats.iter[iTOTAL];
  1120.  if (inc_mode){
  1121.     dc0 = brh->m0.c0 - brh->m1.c0;
  1122.  }else{
  1123.     dc0 = brh->m0.c0;
  1124.  }
  1125.  if (dc0 != 0.){
  1126.     if (brh->n[OUT2].m != 0)
  1127.        recon[brh->n[OUT2].m] += dc0;
  1128.     if (brh->n[OUT1].m != 0)
  1129.        recon[brh->n[OUT1].m] -= dc0;
  1130.  }else if (0){
  1131.     error(bDANGER, "%s: useless load (c0)\n", printlabel(brh,0));
  1132.  }
  1133.  brh->m1 = brh->m0;
  1134. }
  1135. /*--------------------------------------------------------------------------*/
  1136. /* trloadpassive: load passive element to matrix, and offset to right side
  1137.  *     symmetric around diagonal: input == output
  1138.  * if's optimize for speed.  could be done like trloadactive, with OUT for IN.
  1139.  */
  1140. void trloadpassive(brh)
  1141. branch_t *brh;
  1142. {
  1143.  double dc0;
  1144.  double df1;
  1145.  if (brh->loaditer == stats.iter[iTOTAL])
  1146.     error(bDANGER, "%s: double load\n", printlabel(brh,0));
  1147.  brh->loaditer = stats.iter[iTOTAL];
  1148.  if (inc_mode){
  1149.     dc0 = brh->m0.c0 - brh->m1.c0;
  1150.     df1 = brh->m0.f1 - brh->m1.f1;
  1151.  }else{
  1152.     dc0 = brh->m0.c0;
  1153.     df1 = brh->m0.f1;
  1154.  }
  1155.  if (df1 != 0.){
  1156.     if (brh->n[OUT2].m != 0){
  1157.        *AAd(brh->n[OUT2].m,brh->n[OUT2].m) += df1;
  1158.        if (brh->n[OUT1].m != 0){
  1159.       *AAd(brh->n[OUT1].m,brh->n[OUT1].m) += df1;
  1160.       *AA(brh->n[OUT1].m,brh->n[OUT2].m)  -= df1;
  1161.       *AA(brh->n[OUT2].m,brh->n[OUT1].m)  -= df1;
  1162.        }
  1163.     }else if (brh->n[OUT1].m != 0){
  1164.        *AAd(brh->n[OUT1].m,brh->n[OUT1].m) += df1;
  1165.     }
  1166.  }else if (0){
  1167.     error(bDANGER, "%s: useless load (f1)\n", printlabel(brh,0));
  1168.  }
  1169.  if (dc0 != 0.){
  1170.     if (brh->n[OUT2].m != 0)
  1171.        recon[brh->n[OUT2].m] += dc0;
  1172.     if (brh->n[OUT1].m != 0)
  1173.        recon[brh->n[OUT1].m] -= dc0;
  1174.  }else if (0){
  1175.     error(bDANGER, "%s: useless load (c0)\n", printlabel(brh,0));
  1176.  }
  1177.  brh->m1 = brh->m0;
  1178. }
  1179. /*--------------------------------------------------------------------------*/
  1180. /* trloadactive: load active element to matrix, and offset to right side
  1181.  *    asymmetric around diagonal: input != output
  1182.  */
  1183. void trloadactive(brh)
  1184. branch_t *brh;
  1185. {
  1186.  double dc0;
  1187.  double df1;
  1188.  if (brh->loaditer == stats.iter[iTOTAL])
  1189.     error(bDANGER, "%s: double load\n", printlabel(brh,0));
  1190.  brh->loaditer = stats.iter[iTOTAL];
  1191.  if (inc_mode){
  1192.     dc0 = brh->m0.c0 - brh->m1.c0;
  1193.     df1 = brh->m0.f1 - brh->m1.f1;
  1194.  }else{
  1195.     dc0 = brh->m0.c0;
  1196.     df1 = brh->m0.f1;
  1197.  }
  1198.  if (df1 != 0.){
  1199.     *re(brh->n[OUT1].m,brh->n[IN1].m) += df1;
  1200.     *re(brh->n[OUT2].m,brh->n[IN2].m) += df1;
  1201.     *re(brh->n[OUT1].m,brh->n[IN2].m) -= df1;
  1202.     *re(brh->n[OUT2].m,brh->n[IN1].m) -= df1;
  1203.  }else if (0){
  1204.     error(bDANGER, "%s: useless load (f1)\n", printlabel(brh,0));
  1205.  }
  1206.  if (dc0 != 0.){
  1207.     if (brh->n[OUT2].m != 0)
  1208.        recon[brh->n[OUT2].m] += dc0;
  1209.     if (brh->n[OUT1].m != 0)
  1210.        recon[brh->n[OUT1].m] -= dc0;
  1211.  }else if (0){
  1212.     error(bDANGER, "%s: useless load (c0)\n", printlabel(brh,0));
  1213.  }
  1214.  brh->m1 = brh->m0;
  1215. }
  1216. /*--------------------------------------------------------------------------*/
  1217. /* unloadsource: unload source from right side vector
  1218.  */
  1219. void unloadsource(brh)
  1220. branch_t *brh;
  1221. {
  1222.  double dc0;
  1223.  brh->loaditer = 0;
  1224.  brh->m0.c0 = brh->m0.f1 = 0.;
  1225.  if (inc_mode){
  1226.     inc_mode = BAD;
  1227.     dc0 = brh->m0.c0 - brh->m1.c0;
  1228.  }else{
  1229.     dc0 = brh->m0.c0;
  1230.  }
  1231.  if (dc0 != 0.){
  1232.     if (brh->n[OUT2].m != 0)
  1233.        recon[brh->n[OUT2].m] += dc0;
  1234.     if (brh->n[OUT1].m != 0)
  1235.        recon[brh->n[OUT1].m] -= dc0;
  1236.  }
  1237.  brh->m1 = brh->m0;
  1238. }
  1239. /*--------------------------------------------------------------------------*/
  1240. /* unloadpassive: unload passive element and offset
  1241.  *     symmetric around diagonal: input == output
  1242.  * if's optimize for speed.  could be done like unloadactive, with OUT for IN.
  1243.  */
  1244. void unloadpassive(brh)
  1245. branch_t *brh;
  1246. {
  1247.  double dc0;
  1248.  double df1;
  1249.  brh->loaditer = 0;
  1250.  brh->m0.c0 = brh->m0.f1 = 0.;
  1251.  if (inc_mode){
  1252.     inc_mode = BAD;
  1253.     dc0 = brh->m0.c0 - brh->m1.c0;
  1254.     df1 = brh->m0.f1 - brh->m1.f1;
  1255.  }else{
  1256.     dc0 = brh->m0.c0;
  1257.     df1 = brh->m0.f1;
  1258.  }
  1259.  if (df1 != 0.){
  1260.     if (brh->n[OUT2].m != 0){
  1261.        *AAd(brh->n[OUT2].m,brh->n[OUT2].m) += df1;
  1262.        if (brh->n[OUT1].m != 0){
  1263.       *AAd(brh->n[OUT1].m,brh->n[OUT1].m) += df1;
  1264.       *AA(brh->n[OUT1].m,brh->n[OUT2].m)  -= df1;
  1265.       *AA(brh->n[OUT2].m,brh->n[OUT1].m)  -= df1;
  1266.        }
  1267.     }else if (brh->n[OUT1].m != 0){
  1268.        *AAd(brh->n[OUT1].m,brh->n[OUT1].m) += df1;
  1269.     }
  1270.  }
  1271.  if (dc0 != 0.){
  1272.     if (brh->n[OUT2].m != 0)
  1273.        recon[brh->n[OUT2].m] += dc0;
  1274.     if (brh->n[OUT1].m != 0)
  1275.        recon[brh->n[OUT1].m] -= dc0;
  1276.  }
  1277.  brh->m1 = brh->m0;
  1278. }
  1279. /*--------------------------------------------------------------------------*/
  1280. /* unloadactive: unload active element from matrix, and offset from right side
  1281.  *    asymmetric around diagonal: input != output
  1282.  */
  1283. void unloadactive(brh)
  1284. branch_t *brh;
  1285. {
  1286.  double dc0;
  1287.  double df1;
  1288.  brh->loaditer = 0;
  1289.  brh->m0.c0 = brh->m0.f1 = 0.;
  1290.  if (inc_mode){
  1291.     inc_mode = BAD;
  1292.     dc0 = brh->m0.c0 - brh->m1.c0;
  1293.     df1 = brh->m0.f1 - brh->m1.f1;
  1294.  }else{
  1295.     dc0 = brh->m0.c0;
  1296.     df1 = brh->m0.f1;
  1297.  }
  1298.  if (df1 != 0.){
  1299.     *re(brh->n[OUT1].m,brh->n[IN1].m) += df1;
  1300.     *re(brh->n[OUT2].m,brh->n[IN2].m) += df1;
  1301.     *re(brh->n[OUT1].m,brh->n[IN2].m) -= df1;
  1302.     *re(brh->n[OUT2].m,brh->n[IN1].m) -= df1;
  1303.  }
  1304.  if (dc0 != 0.){
  1305.     if (brh->n[OUT2].m != 0)
  1306.        recon[brh->n[OUT2].m] += dc0;
  1307.     if (brh->n[OUT1].m != 0)
  1308.        recon[brh->n[OUT1].m] -= dc0;
  1309.  }
  1310.  brh->m1 = brh->m0;
  1311. }
  1312. /*--------------------------------------------------------------------------*/
  1313. /*--------------------------------------------------------------------------*/
  1314. SHAR_EOF
  1315. fi # end of overwriting check
  1316. if test -f 'src/tr_out.c'
  1317. then
  1318.     echo shar: will not over-write existing file "'src/tr_out.c'"
  1319. else
  1320. cat << \SHAR_EOF > 'src/tr_out.c'
  1321. /* tr_out  01/26/93
  1322.  * Copyright 1983-1992   Albert Davis
  1323.  * tr,dc analysis output functions (and some ac)
  1324.  */
  1325. #include "ecah.h"
  1326. #include "branch.h"
  1327. #include "dc.h"
  1328. #include "io.h"
  1329. #include "mode.h"
  1330. #include "probh.h"
  1331. #include "status.h"
  1332. #include "tr.h"
  1333. #include "worst.h"
  1334. #include "declare.h"
  1335. /*--------------------------------------------------------------------------*/
  1336.     void    dc_out(const dc_t*);
  1337.     void    tr_out(const transient_t*);
  1338.     void    tr_head(double,double,int,char*);
  1339. static    void    tr_print(double);
  1340. /*--------------------------------------------------------------------------*/
  1341. extern const struct ioctrl io;
  1342. extern       struct status stats;
  1343.  
  1344. extern const probe_t *probelist[];
  1345. extern const int sim_mode;
  1346. extern const int worstcase;
  1347. extern const complex_t *fodata;            /* used as a flag */
  1348. /*--------------------------------------------------------------------------*/
  1349. /* dc_out: output the data, "keep" for ac reference
  1350.  * is separate from tr_out only because of old plot functions
  1351.  */
  1352. void dc_out(dc)
  1353. const dc_t *dc;
  1354. {
  1355.  time_start(&(stats.output));
  1356.  if (dc->printnow){
  1357.     trkeep();
  1358.     plotdc(*dc->sweep[0]);
  1359.     tr_print(*dc->sweep[0]);
  1360.     stats.iter[iPRINTSTEP] = 0;
  1361.     stats.control[cSTEPS] = 0;
  1362.  }
  1363.  time_stop(&(stats.output));
  1364. }
  1365. /*--------------------------------------------------------------------------*/
  1366. /* tr_out: output the data, "keep" for ac reference
  1367.  */
  1368. void tr_out(tr)
  1369. const transient_t *tr;
  1370. {
  1371.  time_start(&(stats.output));
  1372.  if (tr->printnow){
  1373.     trkeep();
  1374.     plottr(tr->time);
  1375.     tr_print(tr->time);
  1376.     stats.iter[iPRINTSTEP] = 0;
  1377.     stats.control[cSTEPS] = 0;
  1378.  }
  1379.  time_stop(&(stats.output));
  1380. }
  1381. /*--------------------------------------------------------------------------*/
  1382. /* tr_head: print column headings and draw plot borders
  1383.  */
  1384. void tr_head(start,stop,linear,col1)
  1385. double start;
  1386. double stop;
  1387. int linear;
  1388. char *col1;
  1389. {
  1390.  if (!plopen(sim_mode,start,stop,linear)){
  1391.     if (col1  &&  *col1)
  1392.        mprintf(io.where," %-10s", col1);
  1393.     if (probelist[sim_mode]){
  1394.        int ii;
  1395.        for (ii = 0;  *probelist[sim_mode][ii].what;  ii++){
  1396.       probe_t prb;    /* to hide MSC BUG */
  1397.       prb = probelist[sim_mode][ii];
  1398.       mprintf(io.where," %-10.10s", probename(prb));
  1399.        }
  1400.     }
  1401.     mprintf(io.where,"\n");
  1402.  }
  1403. }
  1404. /*--------------------------------------------------------------------------*/
  1405. /* tr_print: print the list of results (text form) to io.where
  1406.  * The argument is the first column (independent variable, aka "x")
  1407.  */
  1408. static void tr_print(x)
  1409. double x;
  1410. {
  1411.  if (!io.ploton){
  1412.     if (x != NOT_VALID)
  1413.        mprintf(io.where,ftos(x,"           ",5,io.formaat));
  1414.     if (probelist[sim_mode]){
  1415.        int ii;
  1416.        for (ii = 0;  *probelist[sim_mode][ii].what;  ii++){
  1417.       double value;
  1418.       probe_t prb;    /* to hide MSC BUG */
  1419.       prb = probelist[sim_mode][ii];
  1420.       value = trprobe(prb);
  1421.       mprintf(io.where,ftos(value,"           ",5,io.formaat));
  1422.        }
  1423.     }
  1424.     mprintf(io.where,"\n");
  1425.  }
  1426. }
  1427. /*--------------------------------------------------------------------------*/
  1428. /*--------------------------------------------------------------------------*/
  1429. SHAR_EOF
  1430. fi # end of overwriting check
  1431. if test -f 'src/tr_probe.c'
  1432. then
  1433.     echo shar: will not over-write existing file "'src/tr_probe.c'"
  1434. else
  1435. cat << \SHAR_EOF > 'src/tr_probe.c'
  1436. /* tr_probe  01/11/93
  1437.  * Copyright 1983-1992   Albert Davis
  1438.  * transient analysis probe
  1439.  */
  1440. #include "ecah.h"
  1441. #include "branch.h"
  1442. #include "dev.h"
  1443. #include "mode.h"
  1444. #include "nodestat.h"
  1445. #include "probh.h"
  1446. #include "status.h"
  1447. #include "declare.h"
  1448. /*--------------------------------------------------------------------------*/
  1449.     double    trprobe(probe_t);
  1450. static    double    trprobe_node(probe_t);
  1451. /*--------------------------------------------------------------------------*/
  1452. extern const struct status stats;
  1453. extern const struct nodestat *nstat;
  1454. extern const double genout;
  1455. extern const double temp;
  1456. extern const double trtime;
  1457. extern const int *nm;
  1458. extern double *imcon;
  1459. /*--------------------------------------------------------------------------*/
  1460. double trprobe(prb)
  1461. probe_t prb;
  1462. {
  1463.  if (prb.p.node <= stats.total_nodes){
  1464.     return trprobe_node(prb);
  1465.  }else{
  1466.     return trprobe_branch(prb.p.brh,prb.what);
  1467.  }
  1468. }
  1469. /*--------------------------------------------------------------------------*/
  1470. static double trprobe_node(prb)
  1471. probe_t prb;
  1472. {
  1473.  int dummy = 0;
  1474.  
  1475.  setmatch(prb.what,&dummy);
  1476.  if (rematch("ITer")  &&  prb.p.node < iCOUNT){
  1477.     return (double)stats.iter[prb.p.node];
  1478.  }else if (rematch("Control")  &&  prb.p.node < cCOUNT){
  1479.     return (double)stats.control[prb.p.node];
  1480.  }else if (rematch("GEnerator")){
  1481.     return genout;
  1482.  }else if (rematch("Temperature")){
  1483.     return temp + ABS_ZERO;
  1484.  }else if (rematch("TIme")){
  1485.     return trtime;
  1486.  }else if (prb.p.node > stats.total_nodes){
  1487.     return NOT_VALID;
  1488.  }else if (rematch("V")){
  1489.     return imcon[nm[prb.p.node]];
  1490.  }else if (rematch("Logic")){
  1491.     return logicval(nm[prb.p.node]);
  1492.  }else if (rematch("LAstchange")){
  1493.     return nstat[prb.p.node].lastchange;
  1494.  }else if (rematch("FInaltime")){
  1495.     return nstat[prb.p.node].finaltime;
  1496.  }else if (rematch("DIter")){
  1497.     return (double)nstat[prb.p.node].diter;
  1498.  }else if (rematch("AIter")){
  1499.     return (double)nstat[prb.p.node].aiter;
  1500.  }else{ /* bad parameter */
  1501.     return NOT_VALID;
  1502.  }
  1503. }
  1504. /*--------------------------------------------------------------------------*/
  1505. /*--------------------------------------------------------------------------*/
  1506. SHAR_EOF
  1507. fi # end of overwriting check
  1508. if test -f 'src/tr_reviw.c'
  1509. then
  1510.     echo shar: will not over-write existing file "'src/tr_reviw.c'"
  1511. else
  1512. cat << \SHAR_EOF > 'src/tr_reviw.c'
  1513. /* tr_rev  12/31/92
  1514.  * Copyright 1983-1992   Albert Davis
  1515.  * review the solution after solution at a time point
  1516.  * Set up events, evaluate logic inputs, truncation error.
  1517.  * Adjust step size.
  1518.  */
  1519. #include "ecah.h"
  1520. #include "branch.h"
  1521. #include "error.h"
  1522. #include "dev.h"
  1523. #include "mode.h"
  1524. #include "options.h"
  1525. #include "status.h"
  1526. #include "tr.h"
  1527. #include "declare.h"
  1528. /*--------------------------------------------------------------------------*/
  1529.     void    tr_review(transient_t*);
  1530.     double    tr_review_rl(branch_t*);
  1531.     double    tr_review_ll(branch_t**);
  1532. /*--------------------------------------------------------------------------*/
  1533. extern const struct options opt;
  1534. extern         struct status stats;
  1535. extern int sim_phase;    /* simulation phase (init_dc, tran, ...)        */
  1536. /*--------------------------------------------------------------------------*/
  1537. void tr_review(tr)
  1538. transient_t *tr;
  1539. {
  1540.  static double rtime;
  1541.  double tetime;
  1542.  double adt;     /* actual dt */
  1543.  double rdt;    /* review dt */
  1544.  
  1545.  time_start(&(stats.review));
  1546.  ++stats.iter[iTOTAL];
  1547.  if (sim_phase == pINIT_DC){
  1548.     rtime = 0.;
  1549.     adt = rdt = tr->maxdt;
  1550.  }else{
  1551.     rdt = rtime - tr->oldtime;
  1552.     adt = tr->time - tr->oldtime;
  1553.  }
  1554.  if (adt > rdt + opt.mrt)
  1555.     error(bDANGER, "internal error: step control (%g,%g)\n", adt, rdt);
  1556.  
  1557.  if (stats.iter[iSTEP] >= opt.itl[4]){        /* too many iterations */
  1558.     rtime = tr->oldtime + adt / 8.;        /* try again */
  1559.     tr->control = scITER_R;
  1560.  }else if (stats.iter[iSTEP] > opt.itl[3]){    /* too many iterations */
  1561.     rtime = tr->time + adt;            /* no growth */
  1562.     tr->control = scITER_A;
  1563.  }else{                        /* too few iterations */
  1564.     rtime = tr->time + tr->maxdt;        /* use bigger steps */
  1565.     tr->control = scSKIP;
  1566.     if (rtime > tr->time + rdt * 2){
  1567.        rtime = tr->time + rdt * 2;        /* limit to rdt*2 */
  1568.        tr->control = scRDT;
  1569.     }
  1570.     if (rtime > tr->time + adt * 2  &&  rtime > tr->time + rdt){
  1571.        if (rdt > adt * 2)
  1572.       rtime = tr->time + rdt;
  1573.        else                    /* limit to max(rdt,atd*2) */
  1574.       rtime = tr->time + adt * 2;
  1575.        tr->control = scADT;
  1576.     }
  1577.  }
  1578.  tetime = tr_review_rl(firstbranch_dev());    /* trunc error, etc. */
  1579.  time_stop(&(stats.review));
  1580.  if (tetime < rtime){
  1581.     tr->control = scTE;
  1582.     tr->approxtime = tetime;
  1583.  }else{
  1584.     tr->approxtime = rtime;
  1585.  }
  1586. }
  1587. /*--------------------------------------------------------------------------*/
  1588. double tr_review_rl(brh)
  1589. branch_t *brh;
  1590. {
  1591.  double devicetime, worsttime;
  1592.  const branch_t *stop;
  1593.  
  1594.  worsttime = BIGBIG;
  1595.  if (exists(brh)){
  1596.     stop = brh;
  1597.     do {
  1598.        devicetime = tr_review_branch(brh);
  1599.        if (devicetime < worsttime)
  1600.           worsttime = devicetime;
  1601.     } while (brh=nextbranch_dev(brh),  brh != stop);
  1602.  }
  1603.  return worsttime;
  1604. }
  1605. /*--------------------------------------------------------------------------*/
  1606. double tr_review_ll(bl)
  1607. branch_t **bl;
  1608. {
  1609.  double devicetime, worsttime;
  1610.  branch_t *brh;
  1611.  
  1612.  worsttime = BIGBIG;
  1613.  if (bl){
  1614.     for (brh = *bl;  exists(brh);  ++bl){
  1615.        devicetime = tr_review_branch(brh);
  1616.        if (devicetime < worsttime)
  1617.           worsttime = devicetime;
  1618.     }
  1619.  }
  1620.  return worsttime;
  1621. }
  1622. /*--------------------------------------------------------------------------*/
  1623. /*--------------------------------------------------------------------------*/
  1624. SHAR_EOF
  1625. fi # end of overwriting check
  1626. if test -f 'src/tr_setup.c'
  1627. then
  1628.     echo shar: will not over-write existing file "'src/tr_setup.c'"
  1629. else
  1630. cat << \SHAR_EOF > 'src/tr_setup.c'
  1631. /* tr_setup  12/29/92
  1632.  * Copyright 1983-1992   Albert Davis
  1633.  * set up transient and fourier analysis
  1634.  */
  1635. #include "ecah.h"
  1636. #include "argparse.h"
  1637. #include "error.h"
  1638. #include "io.h"
  1639. #include "mode.h"
  1640. #include "options.h"
  1641. #include "tr.h"
  1642. #include "worst.h"
  1643. #include "declare.h"
  1644. /*--------------------------------------------------------------------------*/
  1645.     void    tr_setup(const char*,int*,transient_t*);
  1646.     void    fo_setup(const char*,int*,fourier_t*);
  1647. static    long    to_pow_of_2(double);
  1648. static    void    tr_options(const char*,int*,transient_t*);
  1649. /*--------------------------------------------------------------------------*/
  1650. extern       struct ioctrl io;
  1651. extern const struct options opt;
  1652. extern const double last_time;    /* time to restart at (set by trkeep)        */
  1653. extern int sim_mode;    /* simulation type (AC, DC, ...)            */
  1654. extern int worstcase;    /* worst case, monte carlo mode    (enum type, worst.h)*/
  1655. extern double temp;    /* actual ambient temperature, kelvin            */
  1656. extern int stiff;    /* flag: use stiff integration method            */
  1657. /*--------------------------------------------------------------------------*/
  1658. /* tr_setup: transient analysis: parse command string and set options
  1659.  *     (options set by call to tr_options)
  1660.  */
  1661. void tr_setup(cmd,cnt,tr)
  1662. const char *cmd;
  1663. int  *cnt;
  1664. transient_t *tr;
  1665. {
  1666.  double oldstart, oldstop, oldrange;
  1667.  
  1668.  oldstart = tr->start;
  1669.  oldstop  = tr->stop;
  1670.  oldrange = oldstop - oldstart;
  1671.  
  1672.  tr->cont = YES;
  1673.  if ( isdigit(cmd[*cnt]) || cmd[*cnt]=='.' ){
  1674.     double arg1,arg2,arg3;
  1675.     arg1 = fabs(ctof(cmd,cnt));
  1676.     arg2 = fabs(ctof(cmd,cnt));
  1677.     arg3 = fabs(ctof(cmd,cnt));
  1678.     if (arg3 != 0.){                /* 3 args: all */
  1679.        if (arg1 == 0.  ||  arg1 > arg3){    /* eca (logical) order: */
  1680.       tr->start = arg1;            /* start stop step */
  1681.       tr->stop  = arg2;
  1682.       tr->step  = arg3;
  1683.        }else{                     /* spice (illogical) order */
  1684.       tr->start = arg3;            /* step stop start */
  1685.       tr->stop  = arg2;
  1686.       tr->step  = arg1;
  1687.        }
  1688.     }else if (arg2 != 0.){            /* 2 args */
  1689.        if (arg1 == 0.){                /* 2 args: start, stop */
  1690.       tr->start = arg1;
  1691.       tr->stop  = arg2;
  1692.       /* tr->step unchanged */
  1693.        }else if (arg1 >= arg2){            /* 2 args: stop, step */
  1694.       tr->start = last_time;
  1695.       tr->stop  = arg1;
  1696.       tr->step  = arg2;
  1697.        }else{ /* arg1 < arg2 */            /* 2 args: step, stop */
  1698.       tr->start = arg3; /* 0 */        /* spice order */
  1699.       tr->stop  = arg2;
  1700.       tr->step  = arg1;
  1701.        }
  1702.     }else if (arg1 > last_time){        /* 1 arg: stop */
  1703.        tr->start = last_time;
  1704.        tr->stop  = arg1;
  1705.        /* tr->step unchanged */
  1706.     }else if (arg1 == 0.){            /* 1 arg: start */
  1707.        tr->start = 0.;
  1708.        tr->stop  = oldrange;
  1709.        /* tr->step unchanged */
  1710.     }else{ /* arg1 < last_time, but not 0 */  /* 1 arg: step */
  1711.        tr->start = last_time;
  1712.        tr->stop  = last_time + oldrange;
  1713.        tr->step  = arg1;
  1714.     }
  1715.  }else{ /* no args */
  1716.     tr->start = last_time;
  1717.     tr->stop  = last_time + oldrange;
  1718.     /* tr->step unchanged */
  1719.  }
  1720.  tr_options(cmd,cnt,tr);
  1721.  if  (tr->start < last_time  ||  last_time <= 0.){
  1722.     tr->cont = NO;
  1723.     tr->oldtime = tr->time = 0.;
  1724.  }else{
  1725.     tr->cont = YES;
  1726.     tr->oldtime = tr->time = last_time;
  1727.  }
  1728.  if ( tr->step==0. )
  1729.     error(bERROR, "time step = 0\n");
  1730. }
  1731. /*--------------------------------------------------------------------------*/
  1732. /* fo_setup: fourier analysis: parse command string and set options
  1733.  *     (options set by call to tr_options)
  1734.  */
  1735. void fo_setup(cmd,cnt,fo)
  1736. const char *cmd;
  1737. int  *cnt;
  1738. fourier_t *fo;
  1739. {
  1740.  double oldstop;
  1741.  long points;
  1742.  
  1743.  oldstop = fo->stop;
  1744.  fo->tr.cont = YES;
  1745.  if ( isdigit(cmd[*cnt]) || cmd[*cnt]=='.' ){
  1746.     double arg1,arg2,arg3;
  1747.     arg1 = fabs(ctof(cmd,cnt));
  1748.     arg2 = fabs(ctof(cmd,cnt));
  1749.     arg3 = fabs(ctof(cmd,cnt));
  1750.     if (arg3 != 0.){                /* 3 args: all */
  1751.        fo->start = arg1;
  1752.        fo->stop  = arg2;
  1753.        fo->step  = arg3;
  1754.     }else if (arg2 != 0.){            /* 2 args */
  1755.        if (arg1 >= arg2){            /* 2 args: stop, step */
  1756.       /* fo->start unchanged */
  1757.       fo->stop  = arg1;
  1758.       fo->step  = arg2;
  1759.        }else{ /* arg1 < arg2 */            /* 2 args: start, stop */
  1760.       fo->start = arg1;
  1761.       fo->stop  = arg2;
  1762.       /* fo->step unchanged */
  1763.        }
  1764.     }else if (arg1 > oldstop/8 ){        /* 1 arg: stop */
  1765.        /* fo->start unchanged */
  1766.        fo->stop  = arg1;
  1767.        /* fo->step unchanged */
  1768.     }else if (arg1 == 0.){            /* 1 arg: start */
  1769.     fo->start = 0.;
  1770.     /* fo->stop unchanged */
  1771.     /* fo->step unchanged */
  1772.     }else{ /* arg1 < oldstep/8, but not 0 */  /* 1 arg: step */
  1773.     /* fo->start unchanged */
  1774.     /* fo->stop  unchanged */
  1775.     fo->step  = arg1;
  1776.     }
  1777.  }
  1778.  /* else (no args) : no change */
  1779.  
  1780.  if ( fo->step==0. || fo->stop==0. )
  1781.     error(bERROR, "frequency step = 0\n");
  1782.  tr_options(cmd,cnt,&fo->tr);
  1783.  points = to_pow_of_2(fo->stop*2 / fo->step);
  1784.  fo->tr.time = fo->tr.start = ((fo->tr.cold) ? 0. : last_time);
  1785.  fo->tr.stop = fo->tr.time + 1/fo->step;
  1786.  fo->tr.step = (fo->tr.stop - fo->tr.start) / points;
  1787.  foinit(points);
  1788. }
  1789. /*--------------------------------------------------------------------------*/
  1790. /* to_pow_of_2: round up to nearest power of 2
  1791.  * example: z=92 returns 128
  1792.  */
  1793. static long to_pow_of_2(z)
  1794. double z;
  1795. {
  1796.  long x,y;
  1797.  x = (long)floor(z);
  1798.  for (y=1; x>0; x>>=1 )
  1799.     y<<=1;
  1800.  return y;
  1801. }   
  1802. /*--------------------------------------------------------------------------*/
  1803. /* tr_options: set options common to transient and fourier analysis
  1804.  */
  1805. static void tr_options(cmd,cnt,tr)
  1806. const char *cmd;
  1807. int  *cnt;
  1808. transient_t *tr;
  1809. {
  1810.  io.where |= io.mstdout;
  1811.  temp = opt.tempamb;
  1812.  io.ploton = io.plotset;
  1813.  tr->all = stiff = tr->echo = tr->cold = tr->watch = worstcase = NO;
  1814.  for (;;){
  1815.     if (argparse(cmd,cnt,REPEAT,
  1816.         "ACMAx",    aENUM,        &worstcase,    wMAXAC,
  1817.     "ACMIn",    aENUM,        &worstcase,    wMINAC,
  1818.     "ALl",        aENUM,        &tr->all,    YES,
  1819.     "Ambient",    aODOUBLE,   &temp,    opt.tempamb,
  1820.     "Cold",        aENUM,        &tr->cold,    YES,
  1821.         "DCMAx",    aENUM,        &worstcase,    wMAXDC,
  1822.     "DCMIn",    aENUM,        &worstcase,    wMINDC,
  1823.     "Echo",        aENUM,        &tr->echo,    YES,
  1824.     "LAg",        aENUM,        &worstcase, wLAG,
  1825.     "LEad",        aENUM,        &worstcase, wLEAD,
  1826.     "NOPlot",    aENUM,        &io.ploton, NO,
  1827.     "PLot",        aENUM,        &io.ploton, YES,
  1828.     "RAndom",    aENUM,        &worstcase, wRAND,
  1829.     "Reftemp",    aODOUBLE,   &temp,    opt.tnom,
  1830.     "SEnsitivity",    aENUM,        &worstcase, wSENS,
  1831.     "SKip",        aUINT,        &tr->skip,
  1832.     "STiff",    aENUM,        &stiff,    YES,
  1833.     "Temperature",    aODOUBLE,   &temp,    -ABS_ZERO,
  1834.     "WAtch",    aENUM,        &tr->watch,    YES,
  1835.     ""))
  1836.     ;    
  1837.     else if (outset(cmd,cnt,(char*)NULL,((sim_mode==sTRAN)?"tr ":"fo ")))
  1838.     ;
  1839.     else{
  1840.        syntax(cmd,cnt,bWARNING);
  1841.        break;
  1842.     }
  1843.  }
  1844.  if (tr->skip < 1)
  1845.     tr->skip = 1;
  1846.  tr->maxdt = tr->step / tr->skip;
  1847.  initio(io.where,io.whence);
  1848. }
  1849. /*--------------------------------------------------------------------------*/
  1850. /*--------------------------------------------------------------------------*/
  1851. SHAR_EOF
  1852. fi # end of overwriting check
  1853. if test -f 'src/tr_solve.c'
  1854. then
  1855.     echo shar: will not over-write existing file "'src/tr_solve.c'"
  1856. else
  1857. cat << \SHAR_EOF > 'src/tr_solve.c'
  1858. /* tr_solve  12/29/92
  1859.  * Copyright 1983-1992   Albert Davis
  1860.  * solve one step or a transient analysis
  1861.  */
  1862. #include "ecah.h"
  1863. #include "convstat.h"
  1864. #include "error.h"
  1865. #include "io.h"
  1866. #include "mode.h"
  1867. #include "options.h"
  1868. #include "status.h"
  1869. #include "tr.h"
  1870. #include "declare.h"
  1871. /*--------------------------------------------------------------------------*/
  1872.     int    tr_solve(int);
  1873. /*--------------------------------------------------------------------------*/
  1874. extern       struct ioctrl io;
  1875. extern const struct options opt;
  1876. extern       struct status stats;
  1877.  
  1878. extern int errorcount;    /* count of errors this iteration            */
  1879. extern int inc_mode;    /* make incremental changes                */
  1880. extern int sim_mode;    /* simulation type (AC, DC, ...)            */
  1881. extern int currents_bad;
  1882. /*--------------------------------------------------------------------------*/
  1883. int tr_solve(itl)            /* solve one time step            */
  1884. int itl;
  1885. {
  1886.  int conv;
  1887.  
  1888.  stats.iter[iSTEP] = 0;
  1889.  ++stats.control[cSTEPS];
  1890.  errorcount = 0;
  1891.  
  1892.  for (;;){
  1893.     ++stats.iter[iPRINTSTEP];
  1894.     ++stats.iter[iSTEP];
  1895.     ++stats.iter[sim_mode];
  1896.     ++stats.iter[iTOTAL];
  1897.     tr_clear((inc_mode)?stats.total_nodes:0);
  1898.     conv = tr_fill();
  1899.     if (conv == cGOOD)
  1900.        break;
  1901.     solve();
  1902.     if (stats.iter[iSTEP] >= opt.itl[8]){
  1903.        io.suppresserrors = NO;
  1904.     }
  1905.     if (stats.iter[iSTEP] >= itl){
  1906.        break;
  1907.     }
  1908.     inc_mode = (inc_mode == BAD   ||   stats.iter[iSTEP] == opt.itl[3])
  1909.             ? NO : opt.incmode;
  1910.     currents_bad = NO;
  1911.  }
  1912.  return conv;
  1913. }
  1914. /*--------------------------------------------------------------------------*/
  1915. /*--------------------------------------------------------------------------*/
  1916. SHAR_EOF
  1917. fi # end of overwriting check
  1918. if test -f 'src/tr_sweep.c'
  1919. then
  1920.     echo shar: will not over-write existing file "'src/tr_sweep.c'"
  1921. else
  1922. cat << \SHAR_EOF > 'src/tr_sweep.c'
  1923. /* tr_sweep.c  01/11/93
  1924.  * Copyright 1983-1992   Albert Davis
  1925.  * sweep time and simulate.  output results.
  1926.  * manage event queue
  1927.  */
  1928. #include "ecah.h"
  1929. #include "branch.h"
  1930. #include "convstat.h"
  1931. #include "error.h"
  1932. #include "io.h"
  1933. #include "mode.h"
  1934. #include "options.h"
  1935. #include "status.h"
  1936. #include "tr.h"
  1937. #include "declare.h"
  1938. /*--------------------------------------------------------------------------*/
  1939.     void    tr_sweep(transient_t*);
  1940.     void    new_event(double);
  1941. static    void    tr_first(transient_t*);
  1942. static    int    tr_next(transient_t*);
  1943. /*--------------------------------------------------------------------------*/
  1944. extern       struct ioctrl io;
  1945. extern const struct options opt;
  1946. extern         struct status stats;
  1947. extern const char e_int[];
  1948.  
  1949. extern double trtime;    /* transient analysis time                */
  1950. extern double genout;    /* output of "signal generator" (ckt input)        */
  1951. extern int sim_phase;    /* simulation phase (init_dc, tran, ...)        */
  1952. extern int currents_bad;
  1953.  
  1954. static double events[MAXEVENTCOUNT];    /* the event queue            */
  1955. static int nei;        /* next event index                    */
  1956. static int eeq;        /* end of event queue                    */
  1957. static double nexttick;    /* next pre-scheduled user "event"            */
  1958. /*--------------------------------------------------------------------------*/
  1959. void tr_sweep(tr)
  1960. transient_t *tr;
  1961. {
  1962.  int conv;
  1963.  
  1964.  tr_head(tr->start, tr->stop, YES/*linear*/, "Time");
  1965.  currents_bad = YES;
  1966.  if (tr->cont){            /* use the data from last time        */
  1967.     trestor();
  1968.  }
  1969.  
  1970.  tr_first(tr);
  1971.  sim_phase = pINIT_DC;
  1972.  genout = gen();
  1973.  io.suppresserrors = !tr->watch;
  1974.  conv = tr_solve(opt.itl[1]);
  1975.  if (conv != cGOOD)
  1976.     error(bWARNING, "did not converge\n");
  1977.  tr_review(tr);
  1978.  tr_out(tr);
  1979.  
  1980.  while (tr_next(tr)){
  1981.     sim_phase = pTRAN;
  1982.     genout = gen();
  1983.     io.suppresserrors = !tr->watch;
  1984.     conv = tr_solve(opt.itl[4]);
  1985.     tr_review(tr);
  1986.     if (conv != cGOOD  ||  tr->approxtime <= tr->time)
  1987.        tr->printnow = NO;
  1988.     tr->printnow |= tr->all;
  1989.     tr_out(tr);
  1990.  }
  1991. }
  1992. /*--------------------------------------------------------------------------*/
  1993. void new_event(etime)
  1994. double etime;
  1995. {
  1996.  int ii, jj, kk;
  1997.  
  1998.  /*time_start(&(stats.review));*/
  1999.  for (ii = nei;  ii != eeq;  ii = (ii+1)%MAXEVENTCOUNT){
  2000.     if (etime <= events[ii])
  2001.        break;
  2002.  }
  2003.  for (jj = eeq;  jj != ii;  jj = kk){
  2004.     kk = jj - 1;
  2005.     if (kk < 0)
  2006.        kk += MAXEVENTCOUNT;
  2007.     events[jj] = events[kk];
  2008.  }
  2009.  events[ii] = etime;
  2010.  eeq = (eeq+1)%MAXEVENTCOUNT;
  2011.  if (eeq == nei)
  2012.     error(bERROR, e_int, "event queue overflow");
  2013.  if (opt.foooo == 3){
  2014.     printf("(%d,%d) ", nei, eeq);
  2015.     for (ii = nei;  ii != eeq;  ii = (ii+1)%MAXEVENTCOUNT)
  2016.        printf("%g ", events[ii]);
  2017.     printf("\n");
  2018.  }
  2019.  /*time_stop(&(stats.review));*/
  2020. }
  2021. /*--------------------------------------------------------------------------*/
  2022. static void tr_first(tr)
  2023. transient_t *tr;
  2024. {
  2025.  /* usually, tr->time, tr->oldtime == 0, from tr_setup */
  2026.  time_start(&(stats.review));
  2027.  nexttick = tr->time + tr->step;        /* set next user step */
  2028.  nei = eeq;                    /* empty the queue */
  2029.  trtime = tr->time;
  2030.  tr->printnow = YES;
  2031.  time_stop(&(stats.review));
  2032. }
  2033. /*--------------------------------------------------------------------------*/
  2034. static int tr_next(tr)
  2035. transient_t *tr;
  2036. {
  2037.  double newtime;
  2038.  
  2039.  time_start(&(stats.review));
  2040.  /* tr_review(tr); */                /* trunc error, etc. */
  2041.  
  2042.  newtime = nexttick;                /* user time steps */
  2043.  stats.control[cSTEPCAUSE] = scUSER;
  2044.  
  2045.  if (nei != eeq  &&  events[nei] <= newtime){    /* the event queue */
  2046.     newtime = events[nei];
  2047.     stats.control[cSTEPCAUSE] = scEVENTQ;
  2048.  }
  2049.  
  2050.  if (tr->approxtime < newtime - opt.mrt){
  2051.     if (tr->approxtime < (newtime + tr->time) / 2.){
  2052.        newtime = tr->approxtime;
  2053.        stats.control[cSTEPCAUSE] = tr->control;
  2054.     }else{
  2055.        newtime = (newtime + tr->time) / 2.;
  2056.        stats.control[cSTEPCAUSE] = scHALF;
  2057.     }
  2058.  } 
  2059.  
  2060.  if (nexttick < newtime + opt.mrt){        /* advance user time */
  2061.     tr->printnow = YES;                /* print if user step */
  2062.     nexttick += tr->step;
  2063.  }else{
  2064.     tr->printnow = tr->all;
  2065.  }
  2066.  
  2067.  while (nei != eeq  &&  events[nei] <= newtime){/* advance event queue */
  2068.     nei = (nei+1)%MAXEVENTCOUNT;
  2069.  }
  2070.  
  2071.  if (newtime < tr->time + opt.mrt){
  2072.     if (newtime < tr->time){
  2073.        error(bTRACE, "backwards time step\n");
  2074.        stats.control[cSTEPCAUSE] += scREJECT;
  2075.        if (newtime < (nexttick - tr->step)){
  2076.        error(bTRACE, "user step rejected\n");
  2077.           nexttick -= tr->step;
  2078.        }
  2079.        if (newtime < tr->oldtime){
  2080.           error(bDANGER, "very backward time step\n");
  2081.           newtime = tr->oldtime + opt.mrt;
  2082.        }
  2083.     }else if (newtime == tr->time){
  2084.        error(bWARNING, "zero time step\n");
  2085.        newtime = tr->time + opt.mrt;
  2086.        stats.control[cSTEPCAUSE] += scZERO;
  2087.     }else{
  2088.        error(bWARNING, "time step too small\n");
  2089.        newtime = tr->time + opt.mrt;
  2090.        stats.control[cSTEPCAUSE] += scSMALL;
  2091.     }
  2092.  }
  2093.  
  2094.  if (newtime > tr->time){
  2095.     tr->oldtime = tr->time;
  2096.  }else if (newtime == tr->time){
  2097.     error(bDANGER, e_int, "perpetual zero time step");
  2098.     newtime += opt.mrt;
  2099.     stats.control[cSTEPCAUSE] += scNO_ADVANCE;
  2100.     tr->oldtime = tr->time;        /* if (newtime < tr->time) */
  2101.  }                    /* do not advance oldtime  */
  2102.  
  2103.  trtime = tr->time = newtime;
  2104.  time_stop(&(stats.review));
  2105.  return (tr->time <= tr->stop + opt.mrt);
  2106. }
  2107. /*--------------------------------------------------------------------------*/
  2108. /*--------------------------------------------------------------------------*/
  2109. SHAR_EOF
  2110. fi # end of overwriting check
  2111. #    End of shell archive
  2112. exit 0
  2113.