Index: src/Parser/PcpParser.cpp
===================================================================
--- src/Parser/PcpParser.cpp	(revision 1b2d30e4f20f9cf99960f0856222c0470911beec)
+++ src/Parser/PcpParser.cpp	(revision 650ca8cc5cd79b7bdc9e6dccdeb2e2d235d94356)
@@ -26,5 +26,5 @@
 PcpParser::PcpParser()
 {
-  Parallelization.ProcPEGamma = 1;
+  Parallelization.ProcPEGamma = 8;
   Parallelization.ProcPEPsi = 1;
 
@@ -54,9 +54,9 @@
   LocalizedOrbitals.EpsWannier = 1e-7;
 
-  StepCounts.MaxMinStopStep = 0;
-  StepCounts.InitMaxMinStopStep = 0;
+  StepCounts.MaxMinStopStep = 1;
+  StepCounts.InitMaxMinStopStep = 1;
   StepCounts.OutVisStep = 10;
   StepCounts.OutSrcStep = 5;
-  StepCounts.MaxPsiStep = 3;
+  StepCounts.MaxPsiStep = 0;
   StepCounts.MaxOuterStep = 0;
   StepCounts.MaxMinStep = 100;
@@ -67,5 +67,5 @@
   StepCounts.InitRelEpsTotalEnergy = 1e-05;
   StepCounts.InitRelEpsKineticEnergy = 0.0001;
-  StepCounts.InitMaxMinGapStopStep = 1;
+  StepCounts.InitMaxMinGapStopStep = 0;
 
   PlaneWaveSpecifics.PsiType = 0;
@@ -325,9 +325,21 @@
   class ThermoStatContainer *Thermostats = World::getInstance().getThermostats();
   if (!file->fail()) {
+    // calculate number of Psis
+    vector<atom *> allatoms = World::getInstance().getAllAtoms();
+    CalculateOrbitals(allatoms);
     *file << "# ParallelCarParinello - main configuration file - created with molecuilder" << endl;
     *file << endl;
-    *file << "mainname\t" << Paths.mainname << "\t# programm name (for runtime files)" << endl;
-    *file << "defaultpath\t" << Paths.defaultpath << "\t# where to put files during runtime" << endl;
-    *file << "pseudopotpath\t" << Paths.pseudopotpath << "\t# where to find pseudopotentials" << endl;
+    if (Paths.mainname != NULL)
+      *file << "mainname\t" << Paths.mainname << "\t# programm name (for runtime files)" << endl;
+    else
+      *file << "mainname\tpcp\t# programm name (for runtime files)" << endl;
+    if (Paths.defaultpath != NULL)
+      *file << "defaultpath\t" << Paths.defaultpath << "\t# where to put files during runtime" << endl;
+    else
+      *file << "defaultpath\tnot specified\t# where to put files during runtime" << endl;
+    if (Paths.pseudopotpath != NULL)
+      *file << "pseudopotpath\t" << Paths.pseudopotpath << "\t# where to find pseudopotentials" << endl;
+    else
+      *file << "pseudopotpath\tnot specified\t# where to find pseudopotentials" << endl;
     *file << endl;
     *file << "ProcPEGamma\t" << Parallelization.ProcPEGamma << "\t# for parallel computing: share constants" << endl;
@@ -341,4 +353,6 @@
     *file << "DoFullCurrent\t" << Switches.DoFullCurrent << "\t# Do full perturbation" << endl;
     *file << "DoConstrainedMD\t" << Switches.DoConstrainedMD << "\t# Do perform a constrained (>0, relating to current MD step) instead of unconstrained (0) MD" << endl;
+    ASSERT(Thermostats != NULL, "PcpParser::save() - Thermostats not initialized!");
+    ASSERT(Thermostats->ThermostatNames != NULL, "PcpParser::save() - Thermostats not initialized!");
     *file << "Thermostat\t" << Thermostats->ThermostatNames[Thermostats->Thermostat] << "\t";
     switch(Thermostats->Thermostat) {
@@ -411,14 +425,7 @@
     }
     *file << "PsiType\t\t" << PlaneWaveSpecifics.PsiType << "\t# 0 - doubly occupied, 1 - SpinUp,SpinDown" << endl;
-    // write out both types for easier changing afterwards
-  //  switch (PsiType) {
-  //    case 0:
-        *file << "MaxPsiDouble\t" << PlaneWaveSpecifics.MaxPsiDouble << "\t# here: specifying both maximum number of SpinUp- and -Down-states" << endl;
-  //      break;
-  //    case 1:
-        *file << "PsiMaxNoUp\t" << PlaneWaveSpecifics.PsiMaxNoUp << "\t# here: specifying maximum number of SpinUp-states" << endl;
-        *file << "PsiMaxNoDown\t" << PlaneWaveSpecifics.PsiMaxNoDown << "\t# here: specifying maximum number of SpinDown-states" << endl;
-  //      break;
-  //  }
+    *file << "MaxPsiDouble\t" << PlaneWaveSpecifics.MaxPsiDouble << "\t# here: specifying both maximum number of SpinUp- and -Down-states" << endl;
+    *file << "PsiMaxNoUp\t" << PlaneWaveSpecifics.PsiMaxNoUp << "\t# here: specifying maximum number of SpinUp-states" << endl;
+    *file << "PsiMaxNoDown\t" << PlaneWaveSpecifics.PsiMaxNoDown << "\t# here: specifying maximum number of SpinDown-states" << endl;
     *file << "AddPsis\t\t" << PlaneWaveSpecifics.AddPsis << "\t# Additional unoccupied Psis for bandgap determination" << endl;
     *file << endl;
@@ -427,18 +434,43 @@
     *file << "IsAngstroem\t" << IsAngstroem << "\t# 0 - Bohr, 1 - Angstroem" << endl;
     *file << "RelativeCoord\t" << RelativeCoord << "\t# whether ion coordinates are relative (1) or absolute (0)" << endl;
-    *file << endl;
-    vector<atom *> allatoms = World::getInstance().getAllAtoms();
     map<int, int> ZtoIndexMap;
     OutputElements(file, allatoms, ZtoIndexMap);
     OutputAtoms(file, allatoms, ZtoIndexMap);
-//    result = result && mol->Checkout(file);
-//    if (mol->MDSteps <=1 )
-//      result = result && mol->Output(file);
-//    else
-//      result = result && mol->OutputTrajectories(file);
   } else {
     DoeLog(1) && (eLog()<< Verbose(1) << "Cannot open output file." << endl);
   }
 }
+
+
+/** Counts necessary number of valence electrons and returns number and SpinType.
+ * \param &allatoms all atoms to store away
+ */
+void PcpParser::CalculateOrbitals(vector<atom *> &allatoms)
+{
+  PlaneWaveSpecifics.MaxPsiDouble = PlaneWaveSpecifics.PsiMaxNoDown = PlaneWaveSpecifics.PsiMaxNoUp = PlaneWaveSpecifics.PsiType = 0;
+  for (vector<atom *>::iterator runner = allatoms.begin(); runner != allatoms.end(); ++runner) {
+    PlaneWaveSpecifics.MaxPsiDouble += (*runner)->type->Valence;
+  }
+  cout << PlaneWaveSpecifics.MaxPsiDouble << endl;
+  PlaneWaveSpecifics.PsiMaxNoDown = PlaneWaveSpecifics.MaxPsiDouble/2 + (PlaneWaveSpecifics.MaxPsiDouble % 2);
+  PlaneWaveSpecifics.PsiMaxNoUp = PlaneWaveSpecifics.MaxPsiDouble/2;
+  PlaneWaveSpecifics.MaxPsiDouble /= 2;
+  PlaneWaveSpecifics.PsiType = (PlaneWaveSpecifics.PsiMaxNoDown == PlaneWaveSpecifics.PsiMaxNoUp) ? 0 : 1;
+  if ((PlaneWaveSpecifics.PsiType == 1) && (Parallelization.ProcPEPsi < 2) && ((PlaneWaveSpecifics.PsiMaxNoDown != 1) || (PlaneWaveSpecifics.PsiMaxNoUp != 0))) {
+    Parallelization.ProcPEGamma /= 2;
+    Parallelization.ProcPEPsi *= 2;
+  } else {
+    Parallelization.ProcPEGamma *= Parallelization.ProcPEPsi;
+    Parallelization.ProcPEPsi = 1;
+  }
+  cout << PlaneWaveSpecifics.PsiMaxNoDown << ">" << PlaneWaveSpecifics.PsiMaxNoUp << endl;
+  if (PlaneWaveSpecifics.PsiMaxNoDown > PlaneWaveSpecifics.PsiMaxNoUp) {
+    StepCounts.InitMaxMinStopStep = StepCounts.MaxMinStopStep = PlaneWaveSpecifics.PsiMaxNoDown;
+    cout << PlaneWaveSpecifics.PsiMaxNoDown << " " << StepCounts.InitMaxMinStopStep << endl;
+  } else {
+    StepCounts.InitMaxMinStopStep = StepCounts.MaxMinStopStep = PlaneWaveSpecifics.PsiMaxNoUp;
+    cout << PlaneWaveSpecifics.PsiMaxNoUp << " " << StepCounts.InitMaxMinStopStep << endl;
+  }
+};
 
 /** Prints MaxTypes and list of elements to strea,
@@ -459,4 +491,5 @@
   // print total element count
   *file << "MaxTypes\t" << PresentElements.size() <<  "\t# maximum number of different ion types" << endl;
+  *file << endl;
   // print element list
   *file << "# Ion type data (PP = PseudoPotential, Z = atomic number)" << endl;
@@ -486,5 +519,6 @@
     if (!Inserter.second)
       Inserter.first->second += 1;
-    *file << "Ion_Type" << ZtoIndexMap[(*AtomRunner)->type->Z] << "_" << ZtoCountMap[Inserter.first->second] << "\t"  << fixed << setprecision(9) << showpoint;
+    const int Z = (*AtomRunner)->type->Z;
+    *file << "Ion_Type" << ZtoIndexMap[Z] << "_" << ZtoCountMap[Z] << "\t"  << fixed << setprecision(9) << showpoint;
     *file << (*AtomRunner)->x[0] << "\t" << (*AtomRunner)->x[1] << "\t" << (*AtomRunner)->x[2];
     *file << "\t" << (*AtomRunner)->FixedIon;
