modules/s_sock.cc

   1 /* Copyright (C) 2011-2013 Felix Salfelder, Lars Hedrich
   2  * Author: Felix Salfelder
   3  *
   4  * This file is part of "Gnucap", the Gnu Circuit Analysis Package
   5  *
   6  * This program 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 3, or (at your option)
   9  * any later version.
  10  *
  11  * This program 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.
  15  *
  16  * You should have received a copy of the GNU General Public License
  17  * along with this program; if not, write to the Free Software
  18  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
  19  * 02110-1301, USA.
  20  *------------------------------------------------------------------
  21  * a remote control socket. used by verification tools
  22  */
  23 #include "globals.h"
  24 #include "u_status.h"
  25 #include <unistd.h>
  26 #include "u_prblst.h"
  27 #include "u_cardst.h"
  28 #include "u_nodemap.h"
  29 #include "e_elemnt.h"
  30 #include "e_storag.h"
  31 #include "e_subckt.h"
  32 #include "d_subckt.h"
  33 #include "u_sock.h"
  34 #include "s_ddc.h"
  35 #include "io_error.h"
  36 #include "resolv.h"
  37 #include "s__.h"
  38 #include "io_matrix.h"
  39 #include "m_matrix_extra.h"
  40 #include <iomanip>
  41 #include "d_cap.h"
  42 #include <fcntl.h>
  43 /*--------------------------------------------------------------------------*/
  44 #define userinfo( a,b,c,d,e,f )
  45 /*--------------------------------------------------------------------------*/
  46 using namespace std;
  47 using namespace SOME_CAP_HACK; // FIXME. (maybe use STORAGE device interface?)
  48 /*--------------------------------------------------------------------------*/
  49 namespace {
  50 /*--------------------------------------------------------------------------*/
  51 typedef union {
  52   double double_val;
  53   int32_t   int_val;
  54 } di_union_t;
  55
  56 /*--------------------------------------------------------------------------*/
  57 class SOCK : public DDC_BASE {
  58 public:
  59   explicit SOCK();
  60   ~SOCK();
  61   string status()const;
  62 protected:
  63   void  options(CS&, int x=0);
  64 private:
  65   void  sweep();
  66   void  sweep_recursive(int);
  67   void  first(int);
  68   bool  next(int);
  69   void send_matrix();
  70   void undo_time_step();
  71   void tr_reject();
  72   static unsigned tr_steps_rejected_;
  73 protected:
  74   enum {DCNEST = 4};
  75   int _n_sweeps;
  76 //  typedef void (*p)(double);
  77   //ELEMENT* (_pushel[DCNEST]); /* pointer to thing to sweep, dc command */
  78   CARDSTASH _stash[DCNEST];     /* store std values of elements being swept */
  79   bool _loop[DCNEST];           /* flag: do it again backwards */
  80   bool _reverse_in[DCNEST];     /* flag: sweep backwards, input */
  81   bool _reverse[DCNEST];        /* flag: sweep backwards, working */
  82   bool _cont;                   /* flag: continue from previous run */
  83   TRACE _trace;                 /* enum: show extended diagnostics */
  84   enum {ONE_PT, LIN_STEP, LIN_PTS, TIMES, OCTAVE, DECADE} _stepmode[DCNEST];
  85   static double temp_c_in;      /* ambient temperature, input and sweep variable */
  86   bool _dump_matrix;
  87   unsigned port;
  88   double* U;
  89   double* CU;
  90   double* CUTCU;
  91 public:
  92   void  do_it(CS&, CARD_LIST*);
  93 private:
  94   void  setup(CS&);
  95   void fillnames( const CARD_LIST* scope);
  96   void findcaps( CARD_LIST* scope);
  97   void cap_prepare();
  98   void cap_reset();
  99   vector<string> var_namen_arr;
 100   vector<COMPONENT*> _caplist; // FIXME: use cardlist
 101   CARDSTASH* _capstash;
 102   uint16_t var_namen_total_size;
 103
 104 private: //vera stuff.
 105   void main_loop();
 106   void verainit(unsigned, unsigned, unsigned);
 107   void verakons();
 108   void veraop();
 109   void verainit_reply();
 110   void verakons_reply();
 111   void veraop_reply();
 112   void ev_iter();
 113   void set_param();
 114   unsigned transtep(unsigned init, double dt);
 115   void transtep_reply(unsigned, bool eol=true);
 116   void transtep_gc_reply(unsigned);
 117
 118   char* var_names_buf;
 119
 120   unsigned _verbose;
 121   size_t total;
 122   size_t n_inputs() const{return _input_names.size();}
 123   vector<string> _input_names;
 124   vector<ELEMENT*> _input_devs;
 125   unsigned n_vars;
 126   unsigned n_vars_square;
 127   uint16_t n_eingaenge;
 128   uint16_t length;
 129
 130   SocketStream stream;
 131   //unsigned BUFSIZE;
 132   unsigned n_bytes;
 133
 134   uint16_t _error; // transport between method and method_reply
 135   double _dthack; // same
 136
 137   double *dc_werteA,*dc_loesungA,*kons_loesungA,*kons_residuumA;
 138
 139   int channel;
 140   int frame_number;
 141   Socket* socket;
 142   short unsigned _port_range;
 143   string _port;   // kommt ueber die Uebergabeparamter port
 144                   // globale Variable daher (default: port=1400)
 145   bool _client_mode;
 146   bool _server_mode;
 147   unsigned _bufsize;
 148   bool _bigarg;
 149   string _host;
 150   int reuseaddr;
 151   struct sockaddr_in sin;
 152
 153   double *matrixg, *matrixc,*vectorq;
 154
 155   static const int printlevel=0;
 156
 157 };
 158 /*--------------------------------------------------------------------------*/
 159 unsigned SOCK::tr_steps_rejected_ = 0;
 160 double  SOCK::temp_c_in = 0.;
 161 /*--------------------------------------------------------------------------*/
 162 void SOCK::do_it(CS& Cmd, CARD_LIST* Scope)
 163 {
 164   trace0("SOCK::do_it");
 165   IO::error.detach(stdout);
 166   IO::error.attach(stderr);
 167   _scope = Scope;
 168   _sim->_time0 = 0.;
 169   //_sim->set_command_ddc();
 170   _sim->set_command_dc();
 171   _sim->_phase = p_INIT_DC;
 172   //_sim->_ddc = true;
 173   ::status.ddc.reset().start();
 174   _sim->_temp_c = temp_c_in;
 175   _dump_matrix = 0;
 176   reuseaddr = 0;
 177   _port = "1400";
 178   _port_range = 1;
 179   _client_mode = false;
 180   _host = "localhost";
 181
 182   command_base(Cmd);
 183
 184   //cleanup
 185   cap_reset();
 186 }
 187 /*--------------------------------------------------------------------------*/
 188 SOCK::SOCK() :DDC_BASE(),
 189    _n_sweeps(1),
 190    _cont(false),
 191    _trace(tNONE)
 192 {
 193
 194   for (int ii = 0; ii < DCNEST; ++ii) {
 195     _loop[ii] = false;
 196     _reverse_in[ii] = false;
 197     _reverse[ii] = false;
 198     _step[ii] = 0.;
 199     _linswp[ii] = true;
 200     _zap[ii]=NULL;
 201     _stepmode[ii] = ONE_PT;
 202   }
 203   temp_c_in=OPT::temp_c;
 204   _out=IO::mstdout;
 205 }
 206 /*--------------------------------------------------------------------------*/
 207 void SOCK::setup(CS& Cmd)
 208 {
 209   _sim->_uic = false;
 210   trace0("SOCK::setup");
 211   _cont = false;
 212   _trace = tNONE;
 213   _out = IO::mstdout;
 214   _out.reset(); //BUG// don't know why this is needed */
 215   bool ploton = IO::plotset  &&  plotlist().size() > 0;
 216
 217   options(Cmd);
 218
 219   if (Cmd.more()) {
 220   }else{
 221   }
 222   Cmd.check(bWARNING, "what's this?");
 223
 224   _sim->_uic = true;
 225   CKT_BASE::_sim->init();
 226
 227   IO::plotout = (ploton) ? IO::mstdout : OMSTREAM();
 228   initio(_out);
 229
 230   _error = 0; /* verainit(v_flag, n_inputs, &n_vars, charbuf, &length); */
 231   n_vars = static_cast<uint16_t>( _sim->_total_nodes) ; // _sim->total_nodes doesnt include gnd
 232   // do this also in sweep to initialize the socket streaM
 233   var_namen_arr.resize( n_vars, string("unset"));
 234 //  var_namen_arr[0]="0";
 235   fillnames( &CARD_LIST::card_list );
 236   n_vars_square = n_vars * n_vars;
 237
 238   findcaps(&CARD_LIST::card_list);
 239
 240   // BUG! this kills cap bm's
 241   cap_prepare(); // attach value common if !has_common. stash old common
 242
 243 #ifndef NDEBUG
 244     for (unsigned i=0; i < n_vars; i++)
 245       trace0("name: " + var_namen_arr[i]);
 246 #endif
 247
 248   assert(_n_sweeps > 0);
 249   _sim->_freq = 0;
 250 }
 251 /*--------------------------------------------------------------------------*/
 252 void SOCK::options(CS& Cmd, int Nest)
 253 {
 254   trace0("SOCK::options... ");
 255
 256   _sim->_uic = _loop[Nest] = _reverse_in[Nest] = false;
 257   _port = "1400";
 258   _bufsize = BUFSIZE;
 259   _bigarg = true;
 260   unsigned here = Cmd.cursor();
 261   do{
 262     ONE_OF
 263       || (Cmd.match1("'\"({")   && ((Cmd >> _step_in[Nest]), (_stepmode[Nest] = LIN_STEP)))
 264       || (Cmd.is_float()        && ((Cmd >> _step_in[Nest]), (_stepmode[Nest] = LIN_STEP)))
 265       || (Get(Cmd, "*",           &_step_in[Nest]) && (_stepmode[Nest] = TIMES))
 266       || (Get(Cmd, "+",           &_step_in[Nest]) && (_stepmode[Nest] = LIN_STEP))
 267       || Get(Cmd, "b{ufsize}",    &_bufsize)
 268       || Get(Cmd, "c{ontinue}",   &_cont)
 269       || Get(Cmd, "port" ,        &_port)
 270       || Get(Cmd, "listen{port}", &_port)
 271       || Get(Cmd, "bigarg",       &_bigarg)
 272       || Get(Cmd, "host" ,        &_host)
 273       || Get(Cmd, "tr{s}",        &_do_tran_step)
 274       || Get(Cmd, "dm",           &_dump_matrix)
 275       || Get(Cmd, "client",       &_client_mode)
 276       || Get(Cmd, "server",       &_server_mode)
 277       || Get(Cmd, "dt{emp}",      &temp_c_in,   mOFFSET, OPT::temp_c)
 278       || Get(Cmd, "lo{op}",       &_loop[Nest])
 279       || Get(Cmd, "re{verse}",    &_reverse_in[Nest])
 280       || Get(Cmd, "te{mperature}",&temp_c_in)
 281       || (Cmd.umatch("tr{ace} {=}") &&
 282           (ONE_OF
 283            || Set(Cmd, "n{one}",      &_trace, tNONE)
 284            || Set(Cmd, "o{ff}",       &_trace, tNONE)
 285            || Set(Cmd, "w{arnings}",  &_trace, tUNDER)
 286            || Set(Cmd, "i{terations}",&_trace, tITERATION)
 287            || Set(Cmd, "v{erbose}",   &_trace, tVERBOSE)
 288            || Cmd.warn(bWARNING,
 289                        "need none, off, warnings, iterations, verbose")
 290            )
 291           )
 292       || _out.outset(Cmd);
 293   }while (Cmd.more() && !Cmd.stuck(&here));
 294 }
 295 /*--------------------------------------------------------------------------*/
 296 static void register_status();
 297 /*--------------------------------------------------------------------------*/
 298 void SOCK::sweep()
 299 {
 300   register_status();
 301   // later... FIXME
 302   //
 303   frame_number = 0;
 304   n_bytes = 0;
 305   trace1("SOCK::do_it", _port);
 306
 307   // Calculate buffersize
 308   n_vars = _sim->_total_nodes;
 309   // _sim->total_nodes doesnt include gnd
 310   n_vars_square = n_vars * n_vars;
 311   _bufsize = n_vars_square * 3 * (unsigned)sizeof(double) + 32;
 312   // 3 quadratic matrices * doubles + some safety margin (BUG)
 313
 314   if (_client_mode){
 315     trace1("bufsize Client ", _bufsize);
 316     socket = new ClientSocket(Socket::TCP, _port, _host, _bufsize);
 317     trace1("connected to "+ _host, _port );
 318     stream = *socket;
 319   } else if (_server_mode) {
 320     stringstream p(_port);
 321     uint16_t _port_;
 322     p >> _port_;
 323     trace1("bufsize Server ", _bufsize);
 324     socket = new ServerSocket(Socket::TCP, _port_, 1u, _bufsize);
 325     ServerSocket* sock=prechecked_cast<ServerSocket*>(socket);
 326
 327     trace0("SOCK::do_it waiting");
 328     stream = sock->listen();
 329   } else {
 330     fflush( stdout );
 331     fflush( stdin );
 332     trace0("SOCK::sweep simple i/o");
 333     socket=0;
 334     trace1("bufsize Stdin ", _bufsize);
 335     stream = SocketStream( STDOUT_FILENO, STDIN_FILENO, _bufsize);
 336     stream << "gnucap sock ready";
 337   }
 338
 339   main_loop();
 340   delete socket;
 341   return;
 342 }
 343 /*--------------------------------------------------------------------------*/
 344 string SOCK::status()const
 345 {
 346   // incomplete();
 347   return "sock: tr rejected=" + to_string(tr_steps_rejected_) + "\n";
 348 }
 349 /*--------------------------------------------------------------------------*/
 350 void SOCK::sweep_recursive(int )
 351 {
 352   assert(false);
 353 }
 354 /*--------------------------------------------------------------------------*/
 355 // fetch names from circuit recursively. fill into local vector.
 356 void SOCK::fillnames( const CARD_LIST* scope){
 357   trace0("SOCK::fillnames");
 358
 359   const NODE_MAP * nm = scope->nodes();
 360   for (NODE_MAP::const_iterator i = nm->begin(); i != nm->end(); ++i) {
 361     if (i->first != "0") {
 362       if (const NODE* a= dynamic_cast<const NODE*>(i->second)){
 363         stringstream s;
 364         string myname(a->long_label());
 365         var_namen_arr[a->matrix_number()-1] = myname;
 366         var_namen_total_size = static_cast<uint16_t>( var_namen_total_size + static_cast<uint16_t>(myname.length()) + 1 );
 367       }
 368     }else{
 369     }
 370   }
 371
 372   for (CARD_LIST::const_iterator i = scope->begin(); i != scope->end(); ++i) {
 373     if(const COMPONENT* s = dynamic_cast<const COMPONENT*>(*i)){
 374       if (!s->is_device()){ untested();
 375       }else if ( s->subckt() ) {
 376         fillnames( s->subckt() );
 377       }
 378     }else{ untested();
 379     }
 380   }
 381 }
 382 /*--------------------------------------------------------------------------*/
 383 void SOCK::findcaps( CARD_LIST* scope){
 384   for (CARD_LIST::iterator i = scope->begin(); i != scope->end(); ++i) {
 385     if ( COMPONENT* c = dynamic_cast< COMPONENT*>(*i) ) {
 386       if (c->is_device() && c->has_memory()){
 387         trace1("found cap", c->long_label());
 388         _caplist.push_back( c );
 389       }
 390     }
 391     if (!(*i)->is_device()){ untested();
 392     } else if ( BASE_SUBCKT* s = dynamic_cast< BASE_SUBCKT*>(*i) ) { untested();
 393       trace1("going down", s->long_label());
 394       findcaps( s->subckt() );
 395     }
 396   }
 397 }
 398 /*--------------------------------------------------------------------------*/
 399 static SOCK p2;
 400 static DISPATCHER<CMD>::INSTALL d2(&command_dispatcher, "sock", &p2);
 401 /*--------------------------------------------------------------------------*/
 402 static void register_status()
 403 {
 404   static bool done;
 405
 406   if(!done) {
 407     new DISPATCHER<CKT_BASE>::INSTALL(&status_dispatcher, "sock", &p2);
 408     done = true;
 409   }else{ untested();
 410   }
 411 }
 412 /*--------------------------------------------------------------------------*/
 413 }
 414 /*--------------------------------------------------------------------------*/
 415 static unsigned argc(unsigned opcode)
 416 {
 417   switch(opcode){
 418     case 51: untested();
 419       return 3;
 420     case 52: untested();
 421       return 0;
 422     case 53: untested();
 423       return 0;
 424     case 104: untested();
 425     case 102:
 426       return 1;
 427     default:
 428       return 0;
 429   }
 430 }
 431 /*--------------------------------------------------------------------------*/
 432 void SOCK::main_loop()
 433 {
 434   trace0("SOCK::main_loop");
 435
 436   uint16_t opcode=-1;
 437   uint16_t arg[3];
 438   unsigned char tmp;
 439   arg[0] = -1;
 440   bool init_done=false;
 441
 442   while(true) {
 443     trace0("SOCK::main_loop waiting for opcode");
 444
 445     if(_bigarg){
 446       stream >> tmp >> 7; // sic!
 447       opcode = tmp;
 448       for(unsigned i=0; i<argc(opcode); ++i){
 449         stream >> arg[i] >> 6;
 450       }
 451     }else{
 452       stream >> opcode;
 453       stream >> arg[0];
 454       stream >> arg[1];
 455       stream >> arg[2];
 456     }
 457
 458     _sim->_mode = s_SOCK; // nonsense.
 459                           // use respective built-in mode!
 460                           //
 461     if(_bigarg) {
 462       ::error(bDEBUG, "sock opcode %d\n", opcode);
 463     }else{
 464       ::error(bDEBUG, "sock opcode %d %d %d %d\n", opcode, arg[0], arg[1], arg[2]);
 465     }
 466
 467     double dt;
 468     unsigned status;
 469     switch (opcode) {
 470       case '\0': // 0
 471         return;
 472         break;
 473       case '3': // 51
 474         if(init_done) throw Exception("init twice??");
 475         verainit(arg[0], arg[1], arg[2]);
 476         verainit_reply();
 477         init_done=true;
 478         break;
 479       case '4':  // 52
 480         veraop();
 481         veraop_reply();
 482         break;
 483       case '5':  // 53
 484         verakons();
 485         verakons_reply();
 486         break;
 487       case 101: untested();
 488         ev_iter();
 489         break;
 490       case 102:
 491         stream >> dt;
 492         status = transtep(arg[0], dt);
 493         transtep_reply(status);
 494         break;
 495       case 103: untested();
 496         set_param();
 497         break;
 498       case 104: itested();
 499         stream >> dt;
 500         status = transtep(arg[0], dt);
 501         transtep_gc_reply(status);
 502         break;
 503       default:
 504         ::error(bDANGER, "unknown opcode %i\n", opcode);
 505         throw Exception("unknown opcode");
 506         break;
 507     }
 508     _sim->reset_iteration_counter(iPRINTSTEP); // used by solve_with_homotopy
 509                                                // reset in outdata (not called)
 510   }
 511 }
 512 /*--------------------------------------------------------------------------*/
 513 // very clever way to transfer strings.
 514 static void putstring8(SocketStream* s, const string x)
 515 {
 516   const char* A = x.c_str();
 517
 518   while(*A){
 519     *s<<*A<<*A<<*A<<*A<<*A<<*A<<*A<<*A;
 520     A++;
 521   }
 522 }
 523 /*--------------------------------------------------------------------------*/
 524 void SOCK::verainit(unsigned verbose, unsigned n_in, unsigned length)
 525 { untested();
 526   _verbose = verbose;
 527   char input_namen[length+1];
 528   unsigned here =0;
 529   unsigned n=0;
 530   _input_names.resize(n_in);
 531   _input_devs.resize(n_in);
 532   trace3("verainit: ", verbose,n_inputs(),length);
 533   assert(n_inputs()==n_in);
 534
 535   for (unsigned i=0; i < length; i++) {
 536     stream >> input_namen[i] >> 7;
 537     if(input_namen[i] == '\t'){
 538       input_namen[i] = 0;
 539       trace1("input_namen", string(input_namen+here));
 540       trace5("input_namen",input_namen[i-2],input_namen[i-1], input_namen[i-0],here,i );
 541       _input_names[n] = string(input_namen+here);
 542       here = i+1;
 543
 544       trace1("looking out for", _input_names[n]);
 545       CARD_LIST::fat_iterator ci = findbranch(_input_names[n], &CARD_LIST::card_list);
 546       if (ci.is_end()){
 547         throw Exception("cannot find voltage source \"" + _input_names[n] +"\", giving up");
 548       }
 549       trace0("found input device");
 550       ELEMENT* d = prechecked_cast<ELEMENT*>(*ci);
 551       if (!d){
 552         throw Exception("not something we can use as source, \"" +  _input_names[n] +"\", giving up");
 553       }
 554       _input_devs[n] = d;
 555       assert(_input_devs[n]);
 556
 557       unsigned ii=n;
 558       // see s_cd. FIXME: revert.
 559       _stash[ii] = _input_devs[ii];
 560       _input_devs[ii]->inc_probes();
 561       _input_devs[ii]->set_value(_input_devs[ii]->value(),0);
 562       _input_devs[ii]->set_constant(false);
 563
 564       ++n;
 565     }
 566   }
 567
 568   //trace0("input_namen " + string(input_namen) );
 569   total = (unsigned) (length+4);
 570   assert(stream.bufsize() >= total);
 571
 572   if (!stream.at_end())
 573   {
 574     printf("Error in Verainit! no of bytes received %i <> expected %i\n",
 575         n_bytes, (int)(total*sizeof(di_union_t)));
 576     throw Exception("bloed\n");
 577   }
 578
 579   assert(!var_names_buf);
 580   var_names_buf = (char*) malloc( n_vars * 128 * sizeof(char));
 581   strcpy(var_names_buf,"");
 582   for (unsigned i=0; i < n_vars; i++)
 583   {
 584     trace1("SOCK::verainit ", var_namen_arr[i]);
 585     strcat(var_names_buf, var_namen_arr[i].c_str());
 586     strcat(var_names_buf, "\t");
 587   }
 588   //length = static_cast<uint16_t>( strlen(var_names_buf) );
 589   // userinfo(1,"vera_titan_ak","Variablennamen %s\n",var_names_buf);
 590 }
 591 /*--------------------------------------------------------------------------*/
 592 void SOCK::veraop()
 593 {
 594   trace1("SOCK::veraop",n_vars);
 595   total = n_vars;
 596   _sim->set_command_dc();
 597   assert(stream.bufsize() >= total);
 598
 599   dc_werteA = (double*) malloc(sizeof(double)*n_vars);
 600   trace1("fetching ",n_vars);
 601   assert(_sim->vdc()[0] == 0 );
 602   for (unsigned i=0; i < n_inputs(); i++) {
 603     double d;
 604     stream >> d;
 605     trace3("SOCK::veraop setting input", _input_devs[i]->long_label(), i, d);
 606     _input_devs[i]->set_value(d);
 607   }
 608
 609   _error = 0; /* veraop(sweep_val, x_new, G, C); */
 610
 611   // ================do_dc========
 612   _sim->_uic = false;
 613   _sim->_bypass_ok = false;
 614   _sim->set_inc_mode_bad();
 615   OPT::ITL itl = OPT::DCBIAS;
 616   _sim->_phase = p_INIT_DC;
 617
 618   trace2("SOCK::veraop homotopy", _sim->_phase, _sim->_mode);
 619   if(_dump_matrix) {
 620     _trace=tVERBOSE;
 621   }
 622   CARD_LIST::card_list.tr_begin();  // hier muesste eigentlich eine dc hin.
 623   if(printlist().size()) {
 624     head(0,0," ");
 625   }
 626   bool converged = false;
 627   try{
 628     converged = solve_with_homotopy(itl,_trace);
 629   }catch( Exception e) {
 630     ::error(bDANGER, "hot failed\n");
 631     throw e;
 632   }
 633
 634   if(!converged){ itested();
 635     _error = 1;
 636   }else{
 637
 638     ::status.accept.start();
 639     _sim->set_limit();
 640     trace0("SOCK::veraop, accepting");
 641     CARD_LIST::card_list.tr_accept();
 642     ::status.accept.stop();
 643
 644     _sim->keep_voltages();
 645   }
 646   //========================
 647
 648   // Die Variablenwerte stehen schon durch solve_system in var_werte
 649   // in dc_sysA und muessen noch nach A kopiert werden:
 650   //
 651   //A->var_werte = x_neu;
 652   //free(dc_loesungA);
 653   //
 654   //        A->var_werte = dc_loesungA;
 655   // diff_werte, par_werte muessen auch gesetzt sein!
 656   // in diesem Falle sind sie per default = 0
 657   //
 658   //  solve  somthing
 659   //
 660   //  A->eval_lin_gl(dc_werteA,&matrixg,&matrixc,&vectorq);
 661
 662   trace1("matrix", 1);
 663 }
 664 /*--------------------------------------------------------------------------*/
 665 void SOCK::verakons()
 666 {
 667   //        n_eingaenge == #caps?
 668   _sim->set_command_dc();
 669   trace1("SOCK::kons", _caplist.size());
 670   _sim->_time0 = 0; // just assert?
 671   _sim->_dt0 = 0; // just assert?
 672   advance_time();
 673   _sim->_phase = p_INIT_DC;
 674   total =  n_eingaenge + n_vars + 1;
 675
 676   for (unsigned i=0; i < n_inputs(); i++) {
 677     double d;
 678     stream >> d;
 679     trace3("SOCK::verakons setting input", _input_devs[i]->long_label(), i, d);
 680     _input_devs[i]->set_value(d);
 681   }
 682
 683   assert(_sim->_v0);
 684   for (unsigned i=1; i <= n_vars; i++) {
 685     stream >> _sim->_v0[i];
 686     //    trace2("SOCK::kons start ", i,  _sim->_v0[i] );
 687   }
 688   if (printlist().size()) { untested();
 689     outdata(0.);
 690   }
 691   _sim->keep_voltages(); // v0->vdc
 692
 693   _error = 0; /* verakons(Dwork, x_new, q_dot, G, C); */
 694   //    n_vars = A->n_var;
 695
 696   for( unsigned i = 0; i < _caplist.size(); i++) {
 697     trace1("SOCK::kons",_caplist[i]->long_label());
 698     _caplist[i]->keep_ic(); // latch voltage applied to _v0
 699     _caplist[i]->set_constant(true);
 700     _caplist[i]->q_eval();              // so it will be updated
 701   }
 702   //
 703   //================do_dc========
 704   _sim->_uic = true;
 705   _sim->_bypass_ok = false;
 706   _sim->set_inc_mode_bad();
 707
 708   OPT::ITL itl = OPT::DCBIAS;
 709
 710   if(printlist().size()) { untested();
 711     head(0,0," ");
 712   }
 713   CARD_LIST::card_list.tr_begin();
 714   bool converged = false;
 715   try{
 716     converged = solve(itl,_trace);
 717     //    solve(OPT::TRHIGH,_trace);
 718     //    solve_with_homotopy(itl,_trace);
 719     // homotopy is to much effort calling
 720     // procedures must catch problems with convergence.
 721   }catch( Exception e) {
 722     ::error(bDANGER, "hot failed\n");
 723     throw e;
 724   }
 725   if (!converged) {
 726     ::error(bWARNING, "s_sock::verakons: solve did not converge\n");
 727     converged = solve_with_homotopy(itl,_trace);
 728     if (!converged) {
 729       ::error(bWARNING, "s_sock::verakons: solve did not converge even with homotopy\n");
 730       _error = 1;
 731     }
 732   }
 733   ::status.accept.start();
 734   assert(_sim->_uic);
 735   CARD_LIST::card_list.tr_accept();
 736   ::status.accept.stop();
 737
 738   //  assert(_sim->_mode==s_SOCK);
 739   if (printlist().size()) { untested();
 740     outdata(_sim->_time0);
 741   }
 742   _sim->_mode = s_SOCK; // for now.
 743
 744   _sim->keep_voltages();
 745   _sim->zero_currents();
 746   _sim->set_inc_mode_no();
 747
 748   // vera wants just cap stamps
 749   for (unsigned i = 0; i < _caplist.size(); i++) { itested();
 750     CARD* c = prechecked_cast<CARD*>(_caplist[i]);
 751     assert(c);
 752     trace1("verakons, loading cap", c->long_label());
 753     _sim->_damp = 1.; // need raw stamps
 754     c->tr_load();
 755   }
 756 }
 757 /*--------------------------------------------------------------------------*/
 758 #undef rescale
 759 #define normalize
 760 void SOCK::ev_iter()
 761 {
 762   trace1("SOCK::ev_iter", n_vars);
 763   static double* aug;
 764   double lambda;
 765   if(!aug){ untested(); // quick hack!
 766     assert(_sim->_aa.size()==n_vars);
 767     aug = new double[n_vars*2+1];
 768   }
 769 #ifdef DO_TRACE
 770   for (unsigned i=0; i < 1+2*n_vars; i++){ aug[i]=i; }
 771 #endif
 772   _sim->_aa = -(_sim->_acx.real());
 773   assert(_sim->_aa.size()==n_vars);
 774   double tmp[n_vars+1];
 775   double xiin[n_vars+1];
 776   double norm = 0;
 777   for (unsigned i=0; i < n_vars; i++) {
 778     stream >> tmp[i];
 779     trace1("SOCK::ev_iter got", tmp[i]);
 780     norm+=tmp[i]*tmp[i];
 781   }
 782 #ifdef normalize
 783   norm=1./sqrt(norm);
 784 #endif
 785   for (unsigned i=0; i < n_vars; i++) {
 786 #ifdef normalize
 787     tmp[i]*=norm;
 788 #endif
 789     aug[n_vars*2-i] = tmp[i];
 790     xiin[i]=tmp[i];
 791   }
 792 #ifdef normalize
 793   aug[n_vars] = 0;
 794 #else
 795   aug[n_vars] = norm-1;
 796 #endif
 797   trace1("SOCK filled xi aa", _sim->_aa);
 798   untested();
 799   const BSMATRIX<double> C = _sim->_acx.imag();
 800   C.rmul(aug-1,xiin-1);
 801   stream >> lambda;
 802
 803 #ifdef rescale
 804   for (unsigned i=0; i < n_vars; i++) {aug[i]*=lambda;}
 805 #endif
 806
 807   _sim->_aa+= lambda * _sim->_acx.imag();
 808
 809   trace0("SOCK::ev_iter rmul");
 810   _sim->_aa.rmul(tmp-1,xiin-1);
 811   untested();
 812   tmp[n_vars] = 0;
 813   for (unsigned i=0; i <= n_vars; i++){ untested();
 814     trace2("rhs",i,tmp[i]);
 815   }
 816
 817   _sim->_aa.augment(aug);
 818   assert(_sim->_aa.size()==n_vars+1);
 819   trace1("SOCK augmented aa", _sim->_aa);
 820   _sim->_aa.lu_decomp();
 821   trace1("SOCK decomposed aa", _sim->_aa);
 822   untested();
 823
 824   _sim->_aa.fbsub(tmp-1);
 825   for (unsigned i=0; i <= n_vars; i++){
 826     trace2("res",i,tmp[i]);
 827   }
 828   untested();
 829   trace0("SOCK::ev_iter done");
 830   _sim->_aa.deaugment();
 831   assert(_sim->_aa.size()==n_vars);
 832   stream << ((uint16_t)_error); stream.pad(6);
 833   double resnormsq = 0;
 834   for (unsigned i=0; i < n_vars; i++) { untested();
 835     resnormsq+=tmp[i]*tmp[i];
 836     stream << xiin[i]-tmp[i];
 837   }
 838 #ifdef rescale
 839   resnormsq+=lambda*tmp[n_vars]*lambda*tmp[n_vars];
 840   stream << lambda*(1-tmp[n_vars]);
 841 #else
 842   resnormsq+=tmp[n_vars]*tmp[n_vars];
 843   stream << lambda-tmp[n_vars];
 844 #endif
 845   stream << resnormsq;
 846   stream << SocketStream::eol;
 847   trace0("SOCK::ev_iter sent");
 848 }
 849 /*--------------------------------------------------------------------------*/
 850 void SOCK::set_param()
 851 { untested();
 852   string tmp[2];
 853   stream >> tmp[0];
 854   stream >> tmp[1];
 855   PARAM_LIST* pl = CARD_LIST::card_list.params();
 856   string paramcommand = tmp[0]+"={"+tmp[1]+"}";
 857   trace3("set_param", tmp[0], tmp[1], paramcommand);
 858   CS cmd(CS::_STRING, paramcommand);
 859   pl->parse(cmd);
 860   CARD_LIST::card_list.precalc_first(); // BUG? params don't propagate into sckts otherwise
 861   CARD_LIST::card_list.precalc_last();
 862 }
 863 /*--------------------------------------------------------------------------*/
 864 enum{
 865   sOK=0,          // ok, returned time step indicates next prefered
 866   sNOCONV=1,      // no convergence for second or 3. newton iteration: hope
 867                   // for success with shorter time step
 868   sTRUNC=2,       // Truncation error (time step not accepted (to long!))
 869                   // restart with shorter returned time step
 870   sFIRSTNOCONV=3, // no convergence for first newton iteration after cons
 871                   // time step has no meaning, but reduction of time step should
 872                   // help
 873   sTROUBLE=99     // exception in solve.
 874 };
 875 /*--------------------------------------------------------------------------*/
 876 // init=0 continue with ruse of all old data (do 1 step)
 877 //      1 restore op and do ~5 steps (could be used directly after a not accepted
 878 //        time step)
 879 //      2 do cons_op followed by do ~5 step (depending on Euler,Gear,
 880 //        Trapezoid, options)
 881 //
 882 // return:
 883 // ret (see above)
 884 // dt. positive timestep, use instead in next call.
 885 //
 886 //
 887 unsigned SOCK::transtep(unsigned init, double dt)
 888 {
 889   unsigned ret = sOK;
 890   _error = 0;
 891   trace3("SOCK::transtep", n_vars, init, dt);
 892   _sim->set_command_tran();
 893   _sim->restore_voltages(); //     _vt1[ii] = _v0[ii] = vdc[ii];
 894   _sim->_uic = false;
 895   _sim->_phase = p_RESTORE;
 896
 897   unsigned stepno=-1;
 898   double reftime = _sim->_time0;
 899
 900   if (dt==0.) {
 901     dt = OPT::dtmin;
 902   }else{
 903   }
 904
 905   assert(dt);
 906
 907   if (init==0) { itested();         // continue
 908     stepno = 1;
 909     ::error(bDEBUG, "transient step, at %f, dt %f\n", _sim->_time0, dt);
 910     _sim->_dt0 = dt;
 911     assert(_sim->_time0);
 912     _sim->_time0 += _sim->_dt0;
 913     // CARD_LIST::card_list.tr_restore();
 914   } else if (init==1) { itested(); // restore
 915     assert(_sim->_time0 > 0);
 916     //    CARD_LIST::card_list.tr_restore();
 917     stepno = 1;
 918     _sim->_dt0 = dt/stepno;
 919     ::error(bDEBUG, "restore step, dt %e\n", dt);
 920     _sim->_time0 += _sim->_dt0;
 921   } else if (init==2) {
 922     _sim->_time0 = reftime = 0;
 923     verakons(); // do we need it ?
 924     stepno = 5;
 925     _sim->_dt0 = dt/stepno;
 926     ::error(bDEBUG, "initial step, dt %e\n", dt);
 927     _sim->_time0 += _sim->_dt0;
 928   } else { unreachable();
 929   }
 930
 931
 932   trace2("SOCK::transtep ", _sim->_phase, _sim->_mode);
 933
 934   advance_time(); // sink _time0 into devices
 935   _sim->set_command_tran();
 936
 937   if (printlist().size()) {
 938     outdata(reftime);
 939   }
 940
 941   _sim->_phase = p_TRAN;
 942
 943   //_sim->_bypass_ok = false;
 944   assert(_sim->_dt0 == dt/stepno);
 945   //_sim->_genout = gen();
 946
 947   assert(_sim->analysis_is_tran());
 948
 949
 950   bool tr_converged;
 951   for (unsigned i = stepno; i>0; --i) {
 952     tr_converged = false;
 953     try {
 954       tr_converged = solve(OPT::TRHIGH, _trace);
 955     }catch (Exception e) { incomplete();
 956       ret = sTROUBLE;
 957     }
 958     if (!tr_converged) { itested();
 959       ::error(bDEBUG, "s_sock::transtep did not converge\n");
 960       if (i==stepno) { untested();
 961         ret = sFIRSTNOCONV;
 962         break;
 963       } else { untested();
 964         ret = sNOCONV;
 965         break;
 966       }
 967     }else{
 968     }
 969
 970     if(stepno>1) {
 971       /// uuh ooh HACK.
 972       assert(ret!=sFIRSTNOCONV);
 973       ::status.accept.start(); // accounting
 974       _sim->set_limit();
 975       CARD_LIST::card_list.tr_accept(); // accept solution in all components
 976       ::status.accept.stop();
 977     }
 978     if (printlist().size()) {
 979       outdata(_sim->_time0);
 980     }
 981     if(i>1){
 982       _sim->_time0 += _sim->_dt0;
 983     }
 984   }
 985
 986   TIME_PAIR time_by = CARD_LIST::card_list.tr_review();
 987   double time_by_error_estimate = time_by._error_estimate;
 988   assert(time_by_error_estimate>=0);
 989
 990   // if(!tr_converged?) { ... } else
 991   if (time_by_error_estimate > _sim->_time0) {
 992     dt = time_by_error_estimate - _sim->_time0;
 993     ::status.accept.start();
 994     _sim->set_limit();
 995     CARD_LIST::card_list.tr_accept();
 996     ::status.accept.stop();
 997     _sim->keep_voltages(); //  vdc  = v0
 998     assert(dt>0);
 999   }else{ untested();
1000     double t1 = _sim->_time0 - _sim->_dt0;
1001     dt = time_by_error_estimate - t1;
1002     assert(dt);
1003     ::error(bTRACE,"sock step not accepted at %f, need %f\n", reftime, dt);
1004     ret = sTRUNC;
1005     assert (dt>0); // hopefully.
1006     assert (dt < _sim->_dt0); // hopefully.
1007     dt *= stepno; // stepno is valid for next step.
1008     _sim->_time0 = reftime;
1009     tr_reject();
1010   }
1011
1012 #ifndef NDEBUG
1013   for (unsigned i=1; i <= n_vars; i++) {
1014     if(isnan(_sim->_i[i])||isnan(_sim->_vdcstack.top()[i]) ) {
1015       _error = 1;
1016       assert(!tr_converged);
1017     }
1018   }
1019 #endif
1020
1021   _dthack = dt;
1022   return ret; // FIXME: proper status.
1023 }
1024 /*--------------------------------------------------------------------------*/
1025 void SOCK::tr_reject()
1026 {
1027   ::status.accept.start();
1028   _sim->_acceptq.clear();
1029   ++tr_steps_rejected_;
1030   ::status.accept.stop();
1031 }
1032 /*--------------------------------------------------------------------------*/
1033 void SOCK::verakons_reply()
1034 {
1035   stream << ((uint16_t)_error); stream.pad(6);
1036   for (unsigned i=1; i <= n_vars; i++) {
1037     trace1("SOCK::kons_reply", _sim->vdc()[i]);
1038     stream << _sim->vdc()[i];
1039   }
1040
1041   //
1042   // now i contains the negative sum of the cap value.
1043   // does this break anything?
1044
1045   for (unsigned i=1; i <= n_vars; i++)         /* Q-Punkt == I == RS */
1046   {
1047     trace1("SOCK::verakons_reply  dqdt", _sim->_i[i]);
1048     stream << _sim->_i[i];
1049   }
1050
1051   if(_dump_matrix){
1052     _out << "u,dq=i: \n";
1053     for (unsigned i=0; i <= n_vars; i++)         /* Variablen-Werte */
1054     {
1055       _out << _sim->vdc()[i] << "," << _sim->_i[i] <<  "\n" ;
1056     }
1057   }
1058
1059   trace0("SOCK::verakons_reply  ac_snap");
1060   ac_snapshot(); // FIXME: to verakons()
1061
1062   send_matrix();
1063
1064   assert(stream.bufsize() >= total);
1065   stream << SocketStream::eol;
1066
1067 //  for( unsigned i = 0; i < _caplist.size(); i++)
1068 //  {
1069 //    trace1("unlatching", i);
1070 //
1071 //    dashier umgekehrt.
1072 //    commonstash[i] = _caplist[i]->common();
1073 //    detach(_caplist[i]);
1074 //
1075 //    COMMON_COMPONENT* c = bm_dispatcher.clone("eval_bm_value");
1076 //    COMMON_COMPONENT* dc = c->deflate();
1077 //
1078 //    _caplist[i]->attach_common(dc);
1079 //
1080 //    _caplist[i]->keep_ic(); // latch voltage applied to _v0
1081 //
1082 //
1083 //  }
1084
1085 }
1086
1087 void SOCK::transtep_reply(unsigned ret, bool eol)
1088 {
1089   _error = 0;
1090   stream << _error; stream.pad(6);
1091   stream << (uint64_t)ret;
1092   stream << _dthack;
1093
1094   for (unsigned i=1; i <= n_vars; i++) {
1095     stream << _sim->_vdcstack.top()[i];
1096   }
1097
1098   // OUCH, move to main loop.
1099   if(eol) {
1100     stream << SocketStream::eol;
1101   }
1102 }
1103
1104 void SOCK::transtep_gc_reply(unsigned ret)
1105 { itested();
1106   transtep_reply(ret, false);
1107
1108   for (unsigned i=1; i <= n_vars; i++)
1109   {
1110     trace1("SOCK::transtep_gc_reply  dqdt", _sim->_i[i]);
1111     stream << _sim->_i[i];
1112   }
1113
1114   trace0("SOCK::transtep_gc_reply  ac_snap");
1115   ac_snapshot(); // BUG: calls tr_begin
1116   send_matrix();
1117   assert(stream.bufsize() >= total);
1118   stream << SocketStream::eol;
1119 }
1120 /*--------------------------------------------------------------------------*/
1121 void SOCK::verainit_reply()
1122 {
1123   trace4("SOCK::verainit_reply ", _error, n_vars, length, var_namen_total_size);
1124   stream << _error;   stream.pad(6);
1125   stream << int32_t (n_vars);  stream.pad(4);
1126   stream << int32_t (var_namen_total_size);  stream.pad(4);
1127
1128   assert(stream.tcur() == 24);
1129
1130   for (unsigned i = 0; i < n_vars; i++) {
1131     trace2("SOCK::verainit_reply",  i, var_namen_arr[i]);
1132     putstring8( &stream, var_namen_arr[i]);
1133     stream << '\t'; stream.pad(7);
1134   }
1135
1136 //good idea?? not yet. vera insists on branch node
1137 //
1138 //  for (unsigned i=0; i < n_inputs; i++)
1139 //  {
1140 //    trace1("putting name " +  input_names[i], i);
1141 //    putstring8( &stream, input_names[i]);
1142 //    stream << '\t'; stream.pad(7);
1143 //  }
1144   stream.flush();
1145
1146   trace0("done verainit_reply");
1147 }
1148 /*--------------------------------------------------------------------------*/
1149 void SOCK::veraop_reply()
1150 {
1151   assert(n_vars==_sim->_total_nodes);
1152   trace1("SOCK::veraop_reply ", n_vars);
1153   stream << _error; stream.pad(6);
1154   for (unsigned i=1; i <= n_vars; i++)         /* Variablen-Werte */
1155   {
1156     stream << _sim->vdc()[i];
1157   }
1158   if(_dump_matrix){
1159     _out << "i,u: \n";
1160     for (unsigned i=0; i <= n_vars; i++)         /* Variablen-Werte */
1161     {
1162       _out << _sim->_i[i] << "," <<  _sim->vdc()[i] << "\n" ;
1163     }
1164   }
1165
1166   trace0("SOCK::veraop_reply taking ac snapshot (matrix only?).");
1167   ac_snapshot();
1168   send_matrix();
1169   stream << SocketStream::eol;
1170
1171   trace0("SOCK::veraop_reply sent");
1172
1173   total = 2*n_vars_square+n_vars+1;
1174   assert(stream.bufsize() >= total);
1175   if (printlevel >= 1)
1176   {
1177     userinfo(1,"vera_titan_ak","Sende: Error %i Framenumber %i, Laenge %i\n",
1178         _error,frame_number,total);
1179   }
1180 }
1181 /*--------------------------------------------------------------------------*/
1182 /*--------------------------------------------------------------------------*/
1183 void SOCK::send_matrix()
1184 {
1185   const BSMATRIX<double> G = _sim->_acx.real();
1186   const BSMATRIX<double> C = _sim->_acx.imag();
1187   if(_dump_matrix){
1188      _out << "G\n" << G << "\n";
1189      _out << "C\n" << C << "\n";
1190   }
1191   trace0("SOCK::send_matrix G");
1192   for (unsigned i=1; i <= n_vars; i++){
1193     for (unsigned j=1; j <= n_vars; j++) {
1194       stream << G.s(j,i);
1195     }
1196   }
1197   trace0("SOCK::send_matrix C");
1198   for (unsigned i=1; i <= n_vars; i++){
1199     for (unsigned j=1; j <= n_vars; j++) {
1200       stream << C.s(j,i);
1201     }
1202   }
1203 }
1204 /*--------------------------------------------------------------------------*/
1205 void SOCK::cap_prepare(void)
1206 {
1207   trace1("SOCK::cap_prepare", _caplist.size() );
1208   assert(!_capstash);
1209   _capstash = new CARDSTASH[_caplist.size()];
1210
1211   for (unsigned ii = 0;  ii < _caplist.size();  ++ii) { untested();
1212     _caplist[ii]->inc_probes();                 // we need to keep track of it
1213     _capstash[ii] = _caplist[ii];                       // stash the std value
1214
1215     if(_caplist[ii]->has_common()){ untested();
1216       trace1("SOCK::cap_prepare. have common", _caplist[ii]->value());
1217 //      _caplist[ii]->set_value(_caplist[ii]->value(), 0); // zap out extensions
1218 //      _caplist[ii]->set_constant(false);               // update HACK?
1219 //
1220 //    all devices need individual commons.
1221       _caplist[ii]->attach_common(_caplist[ii]->common()->clone());
1222       assert(_caplist[ii]->has_common());
1223     }else{ untested();
1224       //      _sweepval[ii] = _zap[ii]->set__value();   // point to value to patch
1225       COMMON_COMPONENT* c = bm_dispatcher.clone("eval_bm_value");
1226       double capval = _caplist[ii]->value();
1227       trace2("SOCK::cap_prepare, attaching common", _caplist[ii]->long_label(), capval);
1228       c->set_value(capval);
1229       COMMON_COMPONENT* dc = c->deflate();
1230       assert(dc);
1231       //
1232       _caplist[ii]->set_value(_caplist[ii]->value(), dc);       // zap out extensions
1233       _caplist[ii]->set_constant(false);                // so it will be updated
1234       trace1("SOCK::cap_prepare calling precalc_first", _caplist[ii]->long_label());
1235       _caplist[ii]->precalc_first();
1236       trace0("SOCK::cap_prepare calling precalc_last");
1237       _caplist[ii]->precalc_last();
1238       trace0("SOCK::cap_prepare calling trbegin");
1239       _caplist[ii]->tr_begin();
1240     }
1241   }
1242 }
1243 /*--------------------------------------------------------------------------*/
1244 void SOCK::cap_reset(void)
1245 {
1246   trace0("SOCK::cap_reset");
1247   for (unsigned ii = 0;  ii < _caplist.size();  ++ii) {
1248       _capstash[ii].restore();
1249       _caplist[ii]->dec_probes();
1250 //      _caplist[ii]->precalc_first();
1251 //      _caplist[ii]->precalc_last();
1252   }
1253   delete[] _capstash;
1254   _capstash = NULL;
1255   trace0("SOCK::cap_reset done");
1256 }
1257 /*--------------------------------------------------------------------------*/
1258 SOCK::~SOCK()
1259 {
1260   trace0("SOCK::~SOCK()");
1261 }
1262 /*--------------------------------------------------------------------------*/
1263 // vim:ts=8:sw=2:noet: