}

  }

  {

    xmlChar *v_str = xmlGetProp(cur, BAD_CAST "pin_fixed");

    if (v_str==0) {

      unique_warning("Token \"pin_fixed\" not found in xmlwrite.cpp at or around line no. 1086");

    }

    if (v_str != NULL) {

      pin_fixed = strtobool( (char *)v_str);

      xmlFree(v_str);

    }

  }

  {

    xmlChar *v_str = xmlGetProp(cur, BAD_CAST "at_boundary");

    if (v_str==0) {

      unique_warning("Token \"at_boundary\" not found in xmlwrite.cpp at or around line no. 1097");

    }

    if (v_str != NULL) {

      at_boundary = strtobool( (char *)v_str);

      xmlFree(v_str);

    }

  }

  {

    xmlChar *v_str = xmlGetProp(cur, BAD_CAST "dead");

    if (v_str==0) {

      unique_warning("Token \"dead\" not found in xmlwrite.cpp at or around line no. 1108");

    }

    if (v_str != NULL) {

      dead = strtobool( (char *)v_str);

      xmlFree(v_str);

    }

  }

  {

    xmlChar *v_str = xmlGetProp(cur, BAD_CAST "source");

    if (v_str==0) {

      unique_warning("Token \"source\" not found in xmlwrite.cpp at or around line no. 1119");

    }

    if (v_str != NULL) {

      source = strtobool( (char *)v_str);

      xmlFree(v_str);

    }

  }

  {

    xmlChar *v_str = xmlGetProp(cur, BAD_CAST "boundary");

    if (v_str==0) {

      unique_warning("Token \"boundary\" not found in xmlwrite.cpp at or around line no. 1130");

    }

    if (v_str != NULL) {

      for (int i=0;i<4;i++) {

	if (!xmlStrcmp( v_str, (const xmlChar *)Cell::boundary_type_names[i])) {

	  boundary=(Cell::boundary_type)i;

	  break;

	}

      }

      xmlFree(v_str);

    }

  }

  {

    xmlChar *v_str = xmlGetProp(cur, BAD_CAST "div_counter");

    if (v_str==0) {

      unique_warning("Token \"div_counter\" not found in xmlwrite.cpp at or around line no. 1119");

    }

    if (v_str != NULL) {

      div_counter = atoi( (char *)v_str);

      xmlFree(v_str);

    }

  }

  {

    xmlChar *v_str = xmlGetProp(cur, BAD_CAST "cell_type");

    if (v_str==0) {

      unique_warning("Token \"cell_type\" not found in xmlwrite.cpp at or around line no. 1237");

    }

    if (v_str != NULL) {

      cell_type = atoi( (char *)v_str);

      xmlFree(v_str);

    }

  }

  // Recalculate moments and area

  SetIntegrals();

  return 0;

}

void NodeSet::XMLRead(xmlNode *root, Mesh *m)

{

  xmlNode *n = root->xmlChildrenNode;

  vector<int> tmp_nodes;

  while(n!=NULL) {

    if ((!xmlStrcmp(n->name, (const xmlChar *)"node"))) {

      xmlChar *nc = xmlGetProp(n, BAD_CAST "n");

      if (nc==0) {

	unique_warning("Token \"n\" not found in xmlwrite.cpp at or around line no. 966");

      }

      tmp_nodes.push_back(atoi( (char *)nc));

      xmlFree(nc);

    }

    n = n->next;

  }

  int nnodes = tmp_nodes.size();

  for (int i=0;i<nnodes;i++) {

    AddNode( &(m->getNode(tmp_nodes[i])) );

  }

}

void Wall::XMLAdd(xmlNode *parent) const { 

  // Save the node to a stream so we can reconstruct its state later

  xmlNodePtr xmlwall = xmlNewChild(parent, NULL, BAD_CAST "wall",NULL);

  {

    ostringstream text;

    text << Index();

    xmlNewProp(xmlwall, BAD_CAST "index" , BAD_CAST text.str().c_str() );

  }

  {

    ostringstream text;

    text << c1->Index();

    xmlNewProp(xmlwall, BAD_CAST "c1", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << c2->Index();

    xmlNewProp(xmlwall, BAD_CAST "c2", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << n1->Index();

    xmlNewProp(xmlwall, BAD_CAST "n1", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << n2->Index();

    xmlNewProp(xmlwall, BAD_CAST "n2", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << length;

    xmlNewProp(xmlwall, BAD_CAST "length", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << viz_flux;

    xmlNewProp(xmlwall, BAD_CAST "viz_flux", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << WallTypetoStr(wall_type);

    xmlNewProp(xmlwall, BAD_CAST "wall_type", BAD_CAST text.str().c_str());

  }

  xmlNodePtr tr1_xml = xmlNewChild(xmlwall, NULL, BAD_CAST "transporters1", NULL); 

  if (transporters1) {

    for (int i=0;i<Cell::NChem();i++) {

      xmlNodePtr tr1_val_xml = xmlNewChild(tr1_xml, NULL, BAD_CAST "val", NULL);

      { 

	ostringstream text;

	text << transporters1[i];

	xmlNewProp(tr1_val_xml, BAD_CAST "v", BAD_CAST text.str().c_str());

      }

    }

  }

  if (transporters2) {

    xmlNodePtr tr2_xml = xmlNewChild(xmlwall, NULL, BAD_CAST "transporters2", NULL); 

    for (int i=0;i<Cell::NChem();i++) {

      xmlNodePtr tr2_val_xml = xmlNewChild(tr2_xml, NULL, BAD_CAST "val", NULL);

      { 

	ostringstream text;

	text << transporters2[i];

	xmlNewProp(tr2_val_xml, BAD_CAST "v", BAD_CAST text.str().c_str());

      }

    }

  }

/*  if (apoplast) {

    xmlNodePtr apo_xml = xmlNewChild(xmlwall, NULL, BAD_CAST "apoplast", NULL); 

    for (int i=0;i<Cell::NChem();i++) {

      xmlNodePtr apo_val_xml = xmlNewChild(apo_xml, NULL, BAD_CAST "val", NULL);

      { 

	ostringstream text;

	text << transporters2[i];

	xmlNewProp(apo_val_xml, BAD_CAST "v", BAD_CAST text.str().c_str());

      }

    }

  } */

}

vector<double> XMLIO::XMLReadValArray(xmlNode *cur)

{

  vector<double> result;

  QLocale standardlocale(QLocale::C);

  bool ok;

  xmlNode *valarray_node = cur->xmlChildrenNode;

  while (valarray_node!=NULL) {

    if (!xmlStrcmp(valarray_node->name, (const xmlChar *)"val")) {

      xmlChar *vc = xmlGetProp(valarray_node, BAD_CAST "v");

      if (vc) {

	//result.push_back(strtod( (const char *)vc, 0));

	result.push_back(standardlocale.toDouble((char *)vc, &ok));

	if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(char *)vc);

	xmlFree(vc);

      }

    }

    valarray_node = valarray_node->next;

  }

  return result;

}

void Mesh::XMLSave(const char *docname, xmlNode *options) const

{

  // based on libxml2 example code "tree2.c"

  xmlDocPtr doc = NULL;       /* document pointer */

  xmlNodePtr root_node = NULL;/* node pointers */

  //xmlDtdPtr dtd = NULL;       /* DTD pointer */

  //  LIBXML_TEST_VERSION;

  /* 

   * Creates a new document, a node and set it as a root node

   */

  doc = xmlNewDoc(BAD_CAST "1.0");

  root_node = xmlNewNode(NULL, BAD_CAST "leaf");

  xmlDocSetRootElement(doc, root_node);

  /* 

   * xmlNewProp() creates attributes, which is "attached" to an node.

   * It returns xmlAttrPtr, which isn't used here.

   */

  xmlNewProp(root_node, BAD_CAST "name", BAD_CAST docname);

  time_t t;

  std::time(&t);

  char *tstring = strdup(asctime(localtime(&t))); // but this does

  // replace "end of line character by '\0'

  char *eol=strchr(tstring,'\n');

  if (eol!=NULL)

    *eol='\0';

  xmlNewProp(root_node, BAD_CAST "date", BAD_CAST tstring);

  free(tstring);

  QString simtime = QString("%1").arg(time);

  xmlNewProp(root_node, BAD_CAST "simtime", BAD_CAST simtime.toStdString().c_str()); 

  /*

   * Creates a DTD declaration. Isn't mandatory. 

   */

  par.XMLAdd(root_node);

  xmlNodePtr xmlnodes = xmlNewChild(root_node, NULL, BAD_CAST "nodes",NULL);

  { ostringstream text;

    text << NNodes();

    xmlNewProp(xmlnodes, BAD_CAST "n", BAD_CAST text.str().c_str());

  }

  { ostringstream text;

    text << Node::target_length;

    xmlNewProp(xmlnodes, BAD_CAST "target_length", BAD_CAST text.str().c_str());

  }

  for (vector<Node *>::const_iterator i=nodes.begin(); i!=nodes.end(); i++) {

    (*i)->XMLAdd(xmlnodes) ;

  }

  /* 

   * xmlNewChild() creates a new node, which is "attached" as child node

   * of root_node node. 

   */

  xmlNodePtr xmlcells = xmlNewChild(root_node, NULL, BAD_CAST "cells",NULL);

  {

    ostringstream text;

    text << NCells();

    xmlNewProp(xmlcells, BAD_CAST "n", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << Cell::offset[0];

    xmlNewProp(xmlcells, BAD_CAST "offsetx", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << Cell::offset[1];

    xmlNewProp(xmlcells, BAD_CAST "offsety", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << Cell::factor;

    xmlNewProp(xmlcells, BAD_CAST "magnification", BAD_CAST text.str().c_str());

  }

  {

    ostringstream text;

    text << cells.front()->BaseArea();

    xmlNewProp(xmlcells, BAD_CAST "base_area", BAD_CAST text.str().c_str());

  }

  { 

    ostringstream text;

    text << Cell::NChem();

    xmlNewProp(xmlcells, BAD_CAST "nchem", BAD_CAST text.str().c_str());

  }

  for (vector<Cell *>::const_iterator i=cells.begin(); i!=cells.end(); i++) {

    (*i)->XMLAdd(xmlcells) ;

  }

  boundary_polygon->XMLAdd(xmlcells);

  xmlNodePtr xmlwalls = xmlNewChild(root_node, NULL, BAD_CAST "walls",NULL);

  { 

    ostringstream text;

    text << walls.size();

    xmlNewProp(xmlwalls, BAD_CAST "n", BAD_CAST text.str().c_str());

  }

  for (list<Wall *>::const_iterator i=walls.begin(); i!=walls.end(); i++) {

    (*i)->XMLAdd(xmlwalls) ;

  }

  xmlNodePtr xmlnodesets = xmlNewChild(root_node, NULL, BAD_CAST "nodesets",NULL);

  { 

    ostringstream text;

    text << node_sets.size();

    xmlNewProp(xmlnodesets, BAD_CAST "n", BAD_CAST text.str().c_str());

  }

  for_each( node_sets.begin(), node_sets.end(), bind2nd ( mem_fun( &NodeSet::XMLAdd ), xmlnodesets ) );

  // Add option tree for interactive application

  if (options) {

    xmlAddChild(root_node, options);

  }

  /* 

   * Dumping document to stdio or file

   */

  xmlSaveFormatFileEnc(docname, doc, "UTF-8", 1);

  /*free the document */

  xmlFreeDoc(doc);

  /*

   *Free the global variables that may

   *have been allocated by the parser.

   */

  xmlCleanupParser();

  /*

   * this is to debug memory for regression tests

   */

}

void Mesh::XMLReadSimtime(const xmlNode *a_node)

{

  xmlNode *root_node;

  root_node = (xmlNode *)a_node;

  xmlChar *strsimtime = xmlGetProp(root_node, BAD_CAST "simtime");

  if (strsimtime) {

    QLocale standardlocale(QLocale::C);

    bool ok;

    double simtime=standardlocale.toDouble((char *)strsimtime, &ok);

    if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(char *)strsimtime);

    time = simtime;

#ifdef QDEBUG

    qDebug() << "Simtime = " << strsimtime << endl;

#endif

  } else {

#ifdef QDEBUG

    qDebug() << "No simtime found in file." << endl;

#endif

    time =0;

  }	

}

void Mesh::XMLReadPars(const xmlNode * a_node)

{

  xmlNode *root_node;

  root_node = (xmlNode *)a_node;

  par.XMLRead(root_node);

  Seed(par.rseed);

  MakeDir(par.datadir);

}

void Mesh::XMLReadGeometry(const xmlNode * a_node)

{

  xmlNode *root_node;

  root_node = (xmlNode *)a_node;

  xmlNode *cur;

  QLocale standardlocale(QLocale::C);

  bool ok;

  // allocate Nodes

  cur = root_node->xmlChildrenNode;

  while (cur!=NULL) {

    if ((!xmlStrcmp(cur->name, (const xmlChar *)"nodes"))){

      XMLReadNodes(cur);

    }

    cur=cur->next;

  }

  cur = root_node->xmlChildrenNode;

  while (cur!=NULL) {

    if ((!xmlStrcmp(cur->name, (const xmlChar *)"nodesets"))) {

      XMLReadNodeSets(cur);

    }

    cur=cur->next;

  }

  cur=root_node->xmlChildrenNode;

  while (cur!=NULL) {

    if ((!xmlStrcmp(cur->name, (const xmlChar *)"nodes"))) {

      XMLReadNodeSetsToNodes(cur);

    }

    cur = cur->next;

  }

  // allocate Cells

  cur = root_node;

  // allocate Cells

  cur = cur->xmlChildrenNode;

  while (cur!=NULL) {

    if ((!xmlStrcmp(cur->name, (const xmlChar *)"cells"))){

      xmlChar *offsetxc = xmlGetProp(cur, BAD_CAST "offsetx");

      xmlChar *offsetyc = xmlGetProp(cur, BAD_CAST "offsety");

      double ox=standardlocale.toDouble((const char *)offsetxc, &ok);

      if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)offsetxc);

      double oy=standardlocale.toDouble((const char *)offsetyc, &ok);

      if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)offsetyc);

      Cell::setOffset(ox, oy);

      xmlFree(offsetxc);

      xmlFree(offsetyc);

      xmlChar *magnificationc = xmlGetProp(cur, BAD_CAST "magnification");

      Cell::SetMagnification(standardlocale.toDouble((const char *)magnificationc, &ok));

      if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)magnificationc);

      xmlFree(magnificationc);

      xmlChar *baseareac = xmlGetProp(cur, BAD_CAST "base_area");

      Cell::BaseArea() = standardlocale.toDouble((const char *)baseareac, &ok);

      if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)baseareac);

      xmlFree(baseareac);

      XMLReadCells(cur);

    }

    cur=cur->next;

  }

  // allocate Walls (we need to have read the cells before constructing walls)

  vector <Wall *> tmp_walls;

  cur = root_node->xmlChildrenNode;

  while (cur!=NULL) {

    if ((!xmlStrcmp(cur->name, (const xmlChar *)"walls"))){

      XMLReadWalls(cur, &tmp_walls);

    }

    cur=cur->next;

  }

  // read walls to cells and boundary_polygon

  cur = root_node->xmlChildrenNode;

  while (cur!=NULL) {

    if ((!xmlStrcmp(cur->name, (const xmlChar *)"cells"))) {

      XMLReadWallsToCells(cur, &tmp_walls);

    }

    cur=cur->next;

  }

  boundary_polygon->ConstructNeighborList();

  boundary_polygon->ConstructConnections();

  for (vector<Cell *>::iterator c=cells.begin(); c!=cells.end(); c++) {

    (*c)->ConstructNeighborList();

    (*c)->ConstructConnections();

  }

  shuffled_cells.clear();

  shuffled_cells = cells;

  MyUrand r(shuffled_cells.size());

  random_shuffle(shuffled_cells.begin(),shuffled_cells.end(),r);

}

void Mesh::XMLParseTree(const xmlNode *root_node)

{

  XMLReadSimtime(root_node);

  XMLReadPars(root_node);

  XMLReadGeometry(root_node);

}

void Mesh::XMLReadNodes(xmlNode *root)

{

  QLocale standardlocale(QLocale::C);

  bool ok;

  xmlNode *cur = root;

  cur = cur->xmlChildrenNode;

  for (vector<Node *>::iterator i=nodes.begin(); i!=nodes.end(); i++) {

    delete *i;

  }

  nodes.clear();

  Node::nnodes=0;

  xmlChar *tlc = xmlGetProp(root, BAD_CAST "target_length");

  if (tlc != 0) {

    Node::target_length = standardlocale.toDouble((const char *)tlc, &ok);

    if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)tlc);

    xmlFree(tlc);

  } else {

    // note that libxml2 also defines a token "warning"

    MyWarning::unique_warning("Warning: value found in XML file for Node::target_length.");

  }

  while (cur!=NULL) {

    if ((!xmlStrcmp(cur->name, (const xmlChar *)"node"))){

      xmlChar *xc = xmlGetProp(cur, BAD_CAST "x");

      if (xc==0) {

	unique_warning("Token \"x\" not found in xmlwrite.cpp at or around line no. 722");

      }

      xmlChar *yc = xmlGetProp(cur, BAD_CAST "y");

      if (yc==0) {

	unique_warning("Token \"y\" not found in xmlwrite.cpp at or around line no. 727");

      }

      xmlChar *fixedc = xmlGetProp(cur, BAD_CAST "fixed");

      if (fixedc==0) {

	unique_warning("Token \"fixed\" not found in xmlwrite.cpp at or around line.");

      }

      xmlChar *boundaryc = xmlGetProp(cur, BAD_CAST "boundary");

      if (boundaryc==0) {

	unique_warning("Token \"boundary\" not found in xmlwrite.cpp at or around line.");

      }

      xmlChar *samc = xmlGetProp(cur, BAD_CAST "sam");

      if (samc==0) {

	unique_warning("Token \"sam\" not found in xmlwrite.cpp at or around line.");

      }

      double x = standardlocale.toDouble((const char *)xc, &ok);

      if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)xc);

      double y = standardlocale.toDouble((const char *)yc, &ok);

      if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)yc);

      Node *new_node = new Node(x,y);

      nodes.push_back(new_node);

      new_node->m = this;

      new_node->fixed = strtobool( (char *)fixedc);

      new_node->boundary = strtobool( (char *)boundaryc );

      new_node->sam = strtobool ( (char *)samc);

      new_node->node_set = 0;

      xmlFree(xc);

      xmlFree(yc);

      xmlFree(boundaryc);

      xmlFree(fixedc);

      xmlFree(samc);

    }

    cur=cur->next;

  }

  shuffled_nodes.clear();

  shuffled_nodes = nodes;

  MyUrand r(shuffled_nodes.size());

  random_shuffle(shuffled_nodes.begin(),shuffled_nodes.end(),r);

}

void Mesh::XMLReadWalls(xmlNode *root, vector<Wall *> *tmp_walls)

{

  xmlNode *cur = root;

  cur = cur->xmlChildrenNode;

  for (list<Wall *>::iterator i=walls.begin(); i!=walls.end(); i++) {

    delete *i;

  }

  walls.clear();

  Wall::nwalls = 0;

  tmp_walls->clear();

  QLocale standardlocale(QLocale::C);

  bool ok;

  while (cur!=NULL) {

    vector<int> tmp_nodes;

    while(cur!=NULL) {

      if ((!xmlStrcmp(cur->name, (const xmlChar *)"wall"))) {

	xmlChar *nc = xmlGetProp(cur, BAD_CAST "c1");

	if (nc==0) {

	  unique_warning("Token \"c1\" not found in xmlwrite.cpp at or around line no. 777");

	}

	int c1 = atoi( (char *)nc);

	xmlFree(nc);

	nc = xmlGetProp(cur, BAD_CAST "c2");

	if (	nc==0) {

	  unique_warning("Token \"c2\" not found in xmlwrite.cpp at or around line no. 785");

	}

	int c2 = atoi( (char *)nc);

	xmlFree(nc);

	nc = xmlGetProp(cur, BAD_CAST "n1");

	if (	nc==0) {

	  unique_warning("Token \"n1\" not found in xmlwrite.cpp at or around line no. 793");

	}

	int n1 = atoi( (char *)nc);

	xmlFree(nc);

	nc = xmlGetProp(cur, BAD_CAST "n2");

	if (	nc==0) {

	  unique_warning("Token \"n2\" not found in xmlwrite.cpp at or around line no. 801");

	}

	int n2 = atoi( (char *)nc);

	xmlFree(nc);

	nc = xmlGetProp(cur, BAD_CAST "length");

	if (	nc==0) {

	  unique_warning("Token \"length\" not found in xmlwrite.cpp at or around line no. 809");

	}

	double length = standardlocale.toDouble((const char *)nc, &ok);

	if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)nc);

	xmlFree(nc);

	nc = xmlGetProp(cur, BAD_CAST "viz_flux");

	double viz_flux = 0.0;

	if (nc!=0) {

	  viz_flux = standardlocale.toDouble((const char *)nc, &ok);

	  if (!ok) MyWarning::error("Could not convert \"%s\" to double in XMLRead.",(const char *)nc);

	}

	xmlFree(nc);

	Wall::WallType wall_type;

	{

	  xmlChar *v_str = xmlGetProp(cur, BAD_CAST "wall_type");

	  if (v_str != 0) {

	    if (!xmlStrcmp(v_str, (const xmlChar *)"aux_source")) {

	      wall_type = Wall::AuxSource;

	    } else {

	      if (!xmlStrcmp(v_str, (const xmlChar *)"aux_sink")) {

		wall_type = Wall::AuxSink;

	      } else {

		wall_type = Wall::Normal;

	      }

	    }

	    xmlFree(v_str);

	  } else {

	    wall_type = Wall::Normal;

	  }

	}

	bool dead = false;

	{

	  // Note: property "delete" is used to manually clean up wall lists in XML files

	  // Simply add property 'delete="true"' to the wall and it will be removed from

	  // the mesh. (This saves us from manually reindexing the file). Otherwise do not use it.

	  xmlChar *v_str = xmlGetProp(cur, BAD_CAST "delete");

	  if (v_str != 0) {

	    dead = strtobool( (char *)v_str);

	    xmlFree(v_str);

	  }

	}

	Cell *cc1 = c1 != -1 ? cells[c1] : boundary_polygon;

	Cell *cc2 = c2 != -1 ? cells[c2] : boundary_polygon;

	Wall *w = new Wall( nodes[n1], nodes[n2], cc1, cc2);

	w->length = length;

	w->viz_flux = viz_flux;

	w->wall_type = wall_type;

	w->dead = dead;

	tmp_walls->push_back(w);

	walls.push_back(w);

	xmlNode *w_node = cur->xmlChildrenNode;

	while (w_node!=NULL) {

	  if ((!xmlStrcmp(w_node->name, (const xmlChar *)"transporters1"))) {