     glutMainLoop();
   }
 #endif
-/*
+
   stringstream ss_m;
   ss_m << data_folder << "/" << "timing_" << atoi(argv[1]) << "_" << atoi(argv[2]) << "_" << atoi(argv[3]) << "_" << atoi(argv[4]) << "_" << atof(argv[5]) << ".txt";
   string timing_file_name = ss_m.str();

#ifdef WITH_GLUT
  bool visualize = true;
#endif

  int numElementsPerSide = 4;
  double E = 2e7;


  // if you don't want to visualize, then output the data
  int fileIndex = 0;
  while(sys[workingThread]->time < t_end)
  {
    // Figure out the non-working thread (based on working thread)
    int nonWorkingThread = 1;
    if(workingThread) nonWorkingThread = 0;

    // Output POV-Ray files
    if(sys[workingThread]->getTimeIndex()%outputInterval==0)
    {
      stringstream ss0;
      ss0 << data_folder << "/data0_" << fileIndex << ".dat";
      sys[0]->writeToFile(ss0.str());

      stringstream ss1;
      ss1 << data_folder << "/data1_" << fileIndex << ".dat";
      sys[1]->writeToFile(ss1.str());

      fileIndex++;
    }

    // The working thread will perform the time step while the non-working thread updates the preconditioner
    sys[workingThread]->DoTimeStep();
    if(!sys[nonWorkingThread]->precUpdated) sys[nonWorkingThread]->updatePreconditioner();

    // When the preconditioner is ready, switch the jobs of the systems
    if(sys[workingThread]->timeIndex%200 == 0) {
      cout << "Swap the threads" << endl;
      sys[workingThread]->transferState(sys[nonWorkingThread]);
      workingThread = nonWorkingThread;
    }

    ofile << sys[workingThread]->time                 << ", "
        << sys[workingThread]->stepTime             << ", "
        << sys[workingThread]->stepNewtonIterations << ", "
        << sys[workingThread]->stepKrylovIterations << ", "
        << sys[workingThread]->precUpdated          << " ,     ";
    for (size_t i = 0; i < sys[workingThread]->stepNewtonIterations; ++i)
      ofile << sys[workingThread]->spikeSolveTime[i] << ", " << sys[workingThread]->spikeNumIter[i] << ",     ";
    ofile << endl;
  }

  printf("Total time to simulate: %f [s]\n",sys[workingThread]->timeToSimulate);
  ofile.close();

  return 0;
}

