/*! Type II quasisolitons in Lieb-Liniger: holeIx2

  Purpose:
  Produce the .holeIx2 file
  containing the hole (doubled) quantum numbers of each basis state
  for (multi-)type II hole wavepackets in Lieb-Liniger.

  This executable requires Abacus version 2.

  See README for compilation instructions.

  Copyright © Jean-Sébastien Caux, Anahita Sarvi and Cesare Vianello.

  This program is free software: you can redistribute it and/or modify it under the terms of the GNU Affero General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

  This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Affero General Public License for more details.

  You should have received a copy of the GNU Affero General Public License along with this program. If not, see <https://www.gnu.org/licenses/>.
*/


import std;

import abacus;


int main(int argc, char* argv[])
{
  using namespace std::complex_literals;

  if (argc != 7) {
    std::cout << "Executable holeIx2\n"
	      << " part of the Type II quasisolitons in Lieb-Liniger software suite\n"
	      << " copyright © Jean-Sébastien Caux, Anahita Sarvi and Cesare Vianello.\n";
    std::cout << "\nPurpose:\n"
	      << " Produce the .holeIx2 file\n"
	      << " containing the hole (doubled) quantum numbers of each basis state\n"
	      << " for (multi-)type II hole wavepackets in Lieb-Liniger.\n";

    std::cout << "\nPrerequirements:\n----------------\n"
	      << " - .states file produced by executable rho\n";

    std::cout << "\nUsage:\n------\n";
    std::cout << "holeIx2 <c> <L> <N> <nr holes> <width> <offset>\n\n";

    int warg { 16 }, wtype { 10 }, wcons { 26 };
    std::cout << std::left << std::setw(warg) << "Argument" << std::setw(wtype) << "Type"
	      << std::setw(wcons) << "Constraints" << "Description\n";
    std::cout << std::left << std::setw(warg) << "--------" << std::setw(wtype) << "----"
	      << std::setw(wcons) << "-----------" << "-----------\n";
    std::cout << std::left << std::setw(warg) << "c" << std::setw(wtype) << "Real"
	      << std::setw(wcons) << "> 0"
	      << "Value of the interaction parameter\n";
    std::cout << std::left << std::setw(warg) << "L" << std::setw(wtype) << "Real"
	      << std::setw(wcons) << "> 0"
	      << "System size\n";
    std::cout << std::left << std::setw(warg) << "N" << std::setw(wtype) << "int"
	      << std::setw(wcons) << "> 0"
	      << "Number of particles\n";
    std::cout << std::left << std::setw(warg) << "nr holes" << std::setw(wtype) << "int"
	      << std::setw(wcons) << "1 <= nr holes <= N"
	      << "Number of holes (Type II modes)\n";
    std::cout << std::left << std::setw(warg) << "width" << std::setw(wtype) << "int"
	      << std::setw(wcons) << "nr holes < width <= N"
	      << "Width of the hole window\n";
    std::cout << std::left << std::setw(warg) << "offset" << std::setw(wtype) << "int"
	      << std::setw(wcons) << "0 <= offset <= N-width"
	      << "Offset of the hole window w/r to the right Fermi edge\n";

    return 0;
  }

  std::cout << std::setprecision(std::numeric_limits<Real>::digits10 + 1);

  Real c { std::stold(argv[1]) };
  Real L { std::stold(argv[2]) };
  int N { std::stoi(argv[3]) };
  int nholes { std::stoi(argv[4]) };
  int width { std::stoi(argv[5]) };
  int offset { std::stoi(argv[6]) };


  // Define the input and output files

  std::stringstream filename_base;
  filename_base << "c_" << c << "_N_" << N << "_L_" << L << "_nholes_" << nholes;
  filename_base << "_width_" << width;
  filename_base << "_offset_" << offset;

  std::stringstream states_filename;
  states_filename << filename_base.str() << ".states";
  std::ifstream states_file;
  states_file.open(states_filename.str());
  states_file >> std::setprecision(std::numeric_limits<Real>::digits10 + 1);

  std::stringstream holeIx2_filename;
  holeIx2_filename << filename_base.str() << ".holeIx2";
  std::ofstream holeIx2_file;
  holeIx2_file.open(holeIx2_filename.str(), std::ios::out | std::ios::trunc);
  holeIx2_file << std::setprecision(std::numeric_limits<Real>::digits10 + 1);


  // Process state by state

  std::string tmp_label;
  int tmp_iK;
  Real tmp_E;
  states_file >> tmp_label;
  do {
    states_file >> tmp_iK >> tmp_E;

    ParsedLabel parsed_label(tmp_label);

    holeIx2_file << "\n" << tmp_label;

    // Use the hole index information (hi_) in the parsed plexlabel.
    // Since we assume that labelling is relative to the ground state,
    // the quantum number of the hole is then
    // (hole Ix2) = -(N-1) + 2*(hole index)
    for (int i { 0 }; i < parsed_label.parsed_plexlabel_g_.hi_.size(); ++i)
      holeIx2_file << "\t" << -(N-1) + 2*int(parsed_label.parsed_plexlabel_g_.hi_[i]);

  } while (states_file >> tmp_label);
  states_file.close();

  holeIx2_file.close();

  return 0;
}