process.cpp

#include <string>
#include <cstdio>
#include <iostream>
#include <map>
#include <vector>
#include <boost/algorithm/string.hpp>
#include <array>
using namespace std;
using namespace boost;

struct State
{
    int id;
    int atom_id;
    char atom_name[4];
    int wave_func;
    double j;
    double l;
    double m;
};

// sample state line: state #  1: atom   1 (Sn ), wfc  1 (j=0.5 l=0 m_j=-0.5)
State read_state(const string &stateline)
{

    State state;
    sscanf(stateline.c_str(), "    state #  %d: atom   %d (%s ), wfc  %d (j=%lf l=%lf m_j=%lf)",
           &state.id,
           &state.atom_id,
           state.atom_name,
           &state.wave_func,
           &state.j,
           &state.l,
           &state.m);
    return state;
}

// term example: 0.156*[#429]
struct wf_term
{

    double coeff;
    int state_id;
};

wf_term parse_term(std::string line)
{
    wf_term term;
    sscanf(line.c_str(), "%lf*[#%d]", &term.coeff, &term.state_id);
    return term;
}

struct wave_function
{
    std::vector<wf_term> terms; // the wavefunctin is a sum of terms
};

wave_function read_wavefunction(std::string line)
{

    const string prefix = "psi = ";
    auto match = line.find(prefix);
    if (match != string::npos)
    {

        wave_function wf;
        line = line.substr(match + prefix.size());
        vector<string> terms_str;
        boost::split(terms_str, line, boost::is_any_of("+"));
        for (auto term_expr : terms_str)
        {
            if (term_expr.size() > 0)
            {
                wf.terms.push_back(parse_term(term_expr));
            }
        }
        return wf;
    }
}

typedef vector<pair<double, wave_function>> EigenValueList;
struct KBlock
{

    double k[3];
    EigenValueList eignenvalues;
};

EigenValueList parse_k_block(std::string line)
{
}

int main()
{
    int num_states = 100;
    char line[2000];
    FILE *input_file = fopen("input.txt", "r");

    map<int, State> states;
    for (int i = 0; i < num_states; i++)
    {
        fgets(line, 2000, input_file);
        State st = read_state(line);
        // save in a map for fast lookup
        states[st.id] = st;
    }
    vector<std::pair<array<double, 3>, string>> k_blocks;
    while (fgets(line, 2000, input_file))
    {
        auto l = string(line);
        const string prefix = "k =";
        auto match = l.find(prefix);

        if (match != string::npos)
        {

            double k[3];
            sscanf(l.substr(match + prefix.size()).c_str(), "%lf %lf %lf", k, k + 1, k + 2);
            printf("%lf %lf %lf\n", k[0], k[1], k[2]);
            k_blocks.push_back({{k[0], k[1], k[2]}, ""});
        }
        else if (l.size() > 0 && k_blocks.size())
        {
            auto &block_data = k_blocks.back();
            block_data.second += l;
        }
    }
    for (auto blocks: k_blocks){
        cout << blocks.second << endl;
    }
        //now that I have all the states, read the wavefunctions
        read_wavefunction("  psi = 0.126*[#247]+0.126*[#248]+0.106*[#273]+0.106*[#274]+0.051*[#238]+");
}