text:
//---------------------------------------------------------------------------//
//!
//! \file   tstXsdir.cpp
//! \author Alex Robinson
//! \brief  Xsdir entry unit tests
//!
//---------------------------------------------------------------------------//

// Std Lib Includes
#include <iostream>

// Boost Includes
#include <boost/function.hpp>
#include <boost/bind.hpp>

// Trilinos Includes
#include <Teuchos_UnitTestHarness.hpp>
#include <Teuchos_VerboseObject.hpp>
#include <Teuchos_RCP.hpp>
#include <Teuchos_Array.hpp>

// FRENSIE Includes
#include "DataGen_Xsdir.hpp"
#include "Utility_UnitTestHarnessExtensions.hpp"

//---------------------------------------------------------------------------//
// Testing Variables
//---------------------------------------------------------------------------//
std::string xsdir_file_name;

//---------------------------------------------------------------------------//
// Tests
//---------------------------------------------------------------------------//
// Check that the xsdir file can be exported
TEUCHOS_UNIT_TEST( XsdirEntry, exportInfo )
{
  // Construct the xsdir object from the xsdir file
  DataGen::Xsdir xsdir( xsdir_file_name );

  Teuchos::ParameterList cross_sections( "Cross Sections Directory" );

  xsdir.exportInfo( cross_sections );

  TEST_EQUALITY_CONST( cross_sections.numParams(), 29 );
  TEST_ASSERT( cross_sections.isSublist( "H-1_900.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-1_600.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-1_293.6K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-1_900.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-1_600.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-1_293.6K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-2_900.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-2_600.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-2_293.6K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-2_900.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-2_600.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "H-2_293.6K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-3_900.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-3_600.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-3_293.6K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-3_900.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-3_600.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-3_293.6K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-4_900.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-4_600.0K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-4_293.6K_v8" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-4_900.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-4_600.0K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "He-4_293.6K_v7" ) );
  TEST_ASSERT( cross_sections.isSublist( "be_293.6K_v2" ) );
  TEST_ASSERT( cross_sections.isSublist( "be_400.0K_v2" ) );
  TEST_ASSERT( cross_sections.isSublist( "H" ) );
  TEST_ASSERT( cross_sections.isSublist( "He" ) );
  TEST_ASSERT( cross_sections.isSublist( "Miscellaneous" ) );

  Teuchos::ParameterList& misc_list =
    cross_sections.sublist( "Miscellaneous" );

  TEST_EQUALITY_CONST( misc_list.numParams(), 4 );
  TEST_ASSERT( misc_list.isSublist( "H-2_v7" ) );
  TEST_ASSERT( misc_list.isSublist( "H-2_v2" ) );
  TEST_ASSERT( misc_list.isSublist( "C-12_v7" ) );
  TEST_ASSERT( misc_list.isSublist( "C-12_v2" ) );

  Teuchos::ParameterList& h2_list =
    cross_sections.sublist( "H-2_900.0K_v8" );

  TEST_EQUALITY_CONST( h2_list.numParams(), 21 );

  h2_list = cross_sections.sublist( "H-2_600.0K_v8" );

  TEST_EQUALITY_CONST( h2_list.numParams(), 21 );

  h2_list = cross_sections.sublist( "H-2_293.6K_v8" );

  TEST_EQUALITY_CONST( h2_list.numParams(), 21 );

  h2_list = cross_sections.sublist( "H-2_900.0K_v7" );

  TEST_EQUALITY_CONST( h2_list.numParams(), 21 );

  h2_list = cross_sections.sublist( "H-2_600.0K_v7" );

  TEST_EQUALITY_CONST( h2_list.numParams(), 21 );

  h2_list = cross_sections.sublist( "H-2_293.6K_v7" );

  TEST_EQUALITY_CONST( h2_list.numParams(), 21 );

  Teuchos::ParameterList& other_list =
    cross_sections.sublist( "H-1_900.0K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "H-1_600.0K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "H-1_293.6K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "H-1_900.0K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "H-1_600.0K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "H-1_293.6K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-3_900.0K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-3_600.0K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-3_293.6K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-3_900.0K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-3_600.0K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-3_293.6K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-4_900.0K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-4_600.0K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-4_293.6K_v8" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-4_900.0K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-4_600.0K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "He-4_293.6K_v7" );

  TEST_EQUALITY_CONST( other_list.numParams(), 17 );

  other_list = cross_sections.sublist( "be_293.6K_v2" );

  TEST_EQUALITY_CONST( other_list.numParams(), 5 );

  other_list = cross_sections.sublist( "be_400.0K_v2" );

  TEST_EQUALITY_CONST( other_list.numParams(), 5 );

  other_list = cross_sections.sublist( "H" );

  TEST_EQUALITY_CONST( other_list.numParams(), 10 );

  other_list = cross_sections.sublist( "He" );

  TEST_EQUALITY_CONST( other_list.numParams(), 10 );
}

//---------------------------------------------------------------------------//
// Custom main function
//---------------------------------------------------------------------------//
int main( int argc, char** argv )
{
  Teuchos::CommandLineProcessor& clp = Teuchos::UnitTestRepository::getCLP();

  clp.setOption( "test_xsdir_file",
		 &xsdir_file_name,
		 "Test xsdir file" );

  const Teuchos::RCP<Teuchos::FancyOStream> out =
    Teuchos::VerboseObjectBase::getDefaultOStream();

  Teuchos::CommandLineProcessor::EParseCommandLineReturn parse_return =
    clp.parse(argc,argv);

  if ( parse_return != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL ) {
    *out << "\nEnd Result: TEST FAILED" << std::endl;
    return parse_return;
  }

  // Run the unit tests
  Teuchos::GlobalMPISession mpiSession( &argc, &argv );

  const bool success = Teuchos::UnitTestRepository::runUnitTests( *out );

  if (success)
    *out << "\nEnd Result: TEST PASSED" << std::endl;
  else
    *out << "\nEnd Result: TEST FAILED" << std::endl;

  clp.printFinalTimerSummary(out.ptr());

  return (success ? 0 : 1);
}

//---------------------------------------------------------------------------//
// end tstXsdir.cpp
//---------------------------------------------------------------------------//


text:
/* author: jakob fischer (mail@jakobfischer.eu)
 * description:
 * Tool for solving a ODE for a reaction system given as an jrnf-file. The concentration and time
 * are read and written from / to a comma seperated file. Every row / line represents one time step.
 * the first column contains the time of this step, the second the mean of the quadratic change
 * (<f(x)> ;dx/dt = f(x)). Then follow the concentrations of all species.
 *
 * One version of the program (option "simulate" calculates effective reaction rates directly
 * from formation enthalpies and activation energies given in the network description.
 * (\beta = 1/(k_b T) = 1). This allows higher precision and usage of an explicit solver which
 * is over all faster for small networks (< 50 species).
 *
 * TODO There seems to be a problem with exceptionally short steps (distances in the
 *      <times> vector given to integrate_times) as well as with very wide steps. This
 *      occurs when the "write_log" option lead to output step sizes below 1e-7 or
 *      above 1e7. The solution for now is to first use linear stepping, then
 *      "write_log_abs" and at larger timescale linear stepping again.
 *      Giving a minimal and maximal step size as parameter might be a good
 *      improvement for future versions!
 */

#include <iostream>
#include <fstream>
#include <ctime>
#include <utility>

#include <boost/array.hpp>
#include <boost/numeric/odeint.hpp>

#include <algorithm>
#include "tools/cl_para.h"
#include "net_tools/reaction_network.h"
#include "net_tools/reaction_network_fileop.h"
#include "net_tools/network_tools.h"

using namespace std;
using namespace boost::numeric::odeint;

#ifndef GIT_VERSION
#define GIT_VERSION "no version"
#endif

// version string increment if mayor changes happen
const string odeint_rnet_version_string="1x00x00";

// various ublas types needed for implicit solver
typedef boost::numeric::ublas::vector< double > vector_type;
typedef boost::numeric::ublas::matrix< double > matrix_type;
typedef boost::numeric::ublas::matrix< int > matrix_type_int;

/*
 * Main class that bundles the functionality to load reaction networks and
 * files containing species' concentration. Every object is associated with
 * an network and concentration file at creation. The class offers methods
 * for solving the ODE defined by the network and write the results to the
 * concentration file.
 */

class reaction_network_system {
    double beta;                // beta = 1/(k_b T)
    std::vector<species> sp;    // The reaction
    std::vector<reaction> re;   // network
    // initial time, concentration loaded with concentration file
    double initial_t;
    vector_type initial_con;
    std::string fn_concentration;      // filename

    matrix_type_int N, N_in, N_out;    // stoichiometric matrices
    // List that contains for each column in N_in / N_out (= for each reaction)
    // a vector of those rows (species) that are nonzero.
    std::vector< std::vector<size_t> > N_in_list, N_out_list;
    // Intermediate values / thermodynamic values
    vector_type Ea, mu0, e_bs_in, e_bs_out, e_m_bEa;

public:

    /*
     * Simple calculation of concentration change (rhs / f(x)) from species
     * concentration and reaction constants (<k>, <k_b>).
     */

    struct fast_system {
        reaction_network_system& rns;

        fast_system(reaction_network_system& rns_) : rns(rns_)  {}

        // <x> - state / concentration;  dxdt - change / rhs (output)
        void operator()( const vector_type &x , vector_type &dxdt , double /* t */ ) {

            for(size_t i=0; i<rns.sp.size(); ++i)
                dxdt[i]=0;

            for(size_t i=0; i<rns.re.size(); ++i) {
                double rate_f=rns.re[i].get_k();
	        double rate_b=rns.re[i].get_k_b();

                for(size_t j=0; j<rns.re[i].get_no_educt_s(); ++j)
                    rate_f *= x[rns.re[i].get_educt_id(j)];

                for(size_t j=0; j<rns.re[i].get_no_product_s(); ++j)
                    rate_b *= x[rns.re[i].get_product_id(j)];

                for(size_t j=0; j<rns.re[i].get_no_educt_s(); ++j)
                    dxdt[rns.re[i].get_educt_id(j)] -= (rate_f - rate_b);

                for(size_t j=0; j<rns.re[i].get_no_product_s(); ++j)
                    dxdt[rns.re[i].get_product_id(j)] += (rate_f - rate_b);
             }

            for(size_t i=0; i<rns.sp.size(); ++i)
                if(rns.sp[i].is_constant())
	            dxdt[i]=0;
         }
    };

    /*
     * Calculation of concentration change (rhs) from species concentration
     * and thermodynamic values (energies).
     */

    struct stiff_system {
        reaction_network_system& rns;

        stiff_system(reaction_network_system& rns_) : rns(rns_)  {}

        // Function to calculate reaction rates of state <x> into <rates>.
        void calculate_rates(const vector_type &x, vector_type &rates) {
            vector_type a2(rns.re.size()), a3(rns.re.size());
            for(size_t j=0; j<rns.N.size2(); ++j) {
                a2(j) = 1;
                a3(j) = 1;

                for(size_t k=0; k<x.size(); ++k) {
                    if(rns.N_in(k,j) != 0)
                        a2(j) *= pow(x(k), rns.N_in(k,j));

                    if(rns.N_out(k,j) != 0)
                        a3(j) *= pow(x(k), rns.N_out(k,j));
                }
            }

            rates = element_prod(rns.e_m_bEa, element_prod(a2, rns.e_bs_in) - element_prod(a3, rns.e_bs_out));
        }

        // <x> - state / concentration;  <J> - Jacobi matrix (output);  dfdt - time partial differential (output)
        void operator()( const vector_type &x , vector_type &dxdt , double /* t */ ) {
            vector_type a2(rns.re.size()), a3(rns.re.size());
            for(size_t j=0; j<rns.N.size2(); ++j) {
                a2(j) = 1;
                a3(j) = 1;

                for(size_t l=0; l<rns.N_in_list[j].size(); ++l)
                        a2(j) *= pow(x(rns.N_in_list[j][l]), rns.N_in(rns.N_in_list[j][l],j));

                for(size_t l=0; l<rns.N_out_list[j].size(); ++l)
                        a3(j) *= pow(x(rns.N_out_list[j][l]), rns.N_out(rns.N_out_list[j][l],j));

            }

            dxdt = prod(rns.N, element_prod(rns.e_m_bEa, element_prod(a2, rns.e_bs_in) - element_prod(a3, rns.e_bs_out)));

            for(size_t k=0; k<x.size(); ++k)
                if(rns.sp[k].is_constant())
                    dxdt(k) = 0;
        }
    };

    /*
     * Calculation of Jacobi matrix from species concentration and thermodynamic
     * values (energies).
     */

    struct stiff_system_jacobi {
        const reaction_network_system& rns;

        stiff_system_jacobi(const reaction_network_system& rns_) : rns(rns_)  {}

        // <x> - state / concentration;  <J> - Jacobi matrix (output);  dfdt - time partial differential (output)
        void operator()( const vector_type & x  , matrix_type &J , const double & /* t */ , vector_type &dfdt ) {
            for(size_t i=0; i<x.size(); ++i) {
                for(size_t l=0; l<x.size(); ++l) {
                    double s(0);

                    for(size_t k=0; k<rns.re.size(); ++k) {
                        double p1(1), p2(1);    // p1 = \prod_{j \neq k} {x_j}^N_{jk}^\mathrm{in}
                                                // p2 = \prod_{j \neq k} {x_j}^N_{jk}^\mathrm{out}

                        for(size_t m=0; m<rns.N_in_list[k].size(); ++m) {
                            size_t j=rns.N_in_list[k][m];

                            if(j != l)
                                p1 *= x(j);
                        }

                        for(size_t m=0; m<rns.N_out_list[k].size(); ++m) {
                            size_t j=rns.N_out_list[k][m];

                            if(j != l && rns.N_out(j,k) != 0)
                                p2 *= x(j);
                        }

                        s += rns.N(i,k)*rns.e_m_bEa(k)*(rns.e_bs_in(k)*rns.N_in(l,k)*pow(x(l), rns.N_in(l,k)-1)*p1 -
                                                        rns.e_bs_out(k)*rns.N_out(l,k)*pow(x(l),rns.N_out(l,k)-1)*p2);
                    }

                    J(i,l) = s;
                    if(rns.sp[i].is_constant())
                        J(i,l) = 0;
                }
            }

            dfdt = boost::numeric::ublas::zero_vector<double>(x.size());
        }
    };

    /*
     * Constructor that loads the network file as well as the concentration file
     * (for initial concentration). Also all temporary values that are necessary to
     * solve (calculate rhs + Jacobi matrix) the ode faster are calculated.
     */

    reaction_network_system(const std::string& fn_network, const std::string& fn_concentration_)
        : fn_concentration(fn_concentration_), beta(1) {
        // load network into <re> and <sp>
        read_jrnf_reaction_n(fn_network, sp, re);

	// Remove 1-0 and 0-1 reactions. Such reactions are there to ballance flow through the boundary conditions
        std::remove_if (re.begin(), re.end(),
                        [] (reaction& re) -> bool {
                            return re.get_no_educt() == 1 && re.get_no_product() == 0 ||
                                   re.get_no_educt() == 0 && re.get_no_product() == 1; });

        // Calculate stoichiometric matrices, activation energies and their logarithms...
        N_in = N_out = boost::numeric::ublas::zero_matrix<double>(sp.size(), re.size());
        Ea =  e_bs_in = e_bs_out = e_m_bEa = boost::numeric::ublas::zero_vector<double>(re.size());
        mu0 = boost::numeric::ublas::zero_vector<double>(sp.size());

        // Load species standard gibbs energy / chemical energy into mu0
        for(size_t i=0; i<sp.size(); ++i)
            mu0(i) = sp[i].get_energy();

        // Load / calculate reaction specific quantities like activation Energy
        // <Ea> and joint energy of educts <m0_educts> or products <mu0_products>.
        for(size_t i=0; i<re.size(); ++i) {
            double mu0_educts(0), mu0_products(0);

            for(size_t j=0; j<re[i].get_no_educt_s(); ++j) {
                size_t id(re[i].get_educt_id(j)),  mul(re[i].get_educt_mul(j));
                N_in(id,i) += mul;
                mu0_educts += sp[id].get_energy()*mul;
            }

            for(size_t j=0; j<re[i].get_no_product_s(); ++j) {
                size_t id(re[i].get_product_id(j)),  mul(re[i].get_product_mul(j));
                N_out(id,i) += mul;
                mu0_products += sp[id].get_energy()*mul;
            }

            Ea(i) = re[i].get_activation()+max(mu0_educts,mu0_products);
            e_bs_in(i) = exp(beta*mu0_educts);
            e_bs_out(i) = exp(beta*mu0_products);
            e_m_bEa(i) = exp(-beta*Ea(i));
        }

        N = N_out-N_in;

        // Build lists (for speedup)
        for(size_t i=0; i<re.size(); ++i) {
            N_in_list.push_back(std::vector<size_t>());
            N_out_list.push_back(std::vector<size_t>());

            for(size_t j=0; j<sp.size(); ++j) {
                if(N_in(j,i) != 0)
                    N_in_list[i].push_back(j);

                if(N_out(j,i) != 0)
                    N_out_list[i].push_back(j);
            }
        }

        // Read last line of concentration file and write initial concentration to initial_con
        // and initial time to initial_t. When done, open the same file for appending...
        initial_con = boost::numeric::ublas::zero_vector<double>(sp.size());

        std::ifstream  data(fn_concentration.c_str());

        if(!data.good()) {
           std::cout << "Could not open concentration file: " << fn_concentration << std::endl;
           return;
        }

        std::string line;
        std::getline(data,line);       // dont want the header
        double last_msd=0;

        // Get line of actual data
        while(!std::getline(data,line).eof()) {
            std::stringstream ls(line);
            std::string cell;

            // read in value by value and count its position to place it in the
            // right variables
            size_t cnt=0;
            while(std::getline(ls,cell,',')) {
                std::stringstream in(cell);

                if(cnt == 0)
                    in >> initial_t ;
                else if(cnt == 1)
                    in >> last_msd;
                else if(cnt <= sp.size()+1)
                    in >> initial_con(cnt-2);
                else
                    std::cout << "Error at reading csv / concentration file!" << std::endl;

                ++cnt;
            }
        }

        std::cout << "Simulating file: " << fn_concentration << std::endl;
        std::cout << "Loaded concentration file with starting time " << initial_t << std::endl;
        data.close();
    }

    /*
     * To check some differential equations or the algorithm it might be usefull
     * to calculate the right hand side of the differential equation for a given
     * concentration vector.
     */

    void print_rhs() {
        vector_type x(initial_con);
        vector_type dxdt=boost::numeric::ublas::zero_vector<double>(sp.size());

        // calculate rhs - operator saves it do dxdt
        stiff_system(*this)(x, dxdt, 0);

        cout << "righthandside:" << endl;
        for(size_t i=0; i<dxdt.size(); ++i)
            cout << "/  " << dxdt(i) << "  /";
        cout << endl;
    }

    /*
     * Method for calling the system that uses energetics for calculating effective
     * rates and allows usage of a stiff solver - implicit stepper. System (f(x))
     * is defined in "stiff_system" (rhs) and "stiff_system_jacobi (Jacobi matrix).
     *
     * Tmax           - time up to that ode is solved
     * deltaT         - initial step size for addaptive stepper (does not relate to time output is written)
     * solve_implicit - do use implicit solver?
     * wint           - number of time output is written to output file between time of last
     *                  entry in concentration file and Tmax
     * write_log = 0  - linearly spaced times for output
     *           = 1  - logarithmically spaced (relative -> log(t-t0) is equidist.)
     *           = 2  - logarithmically spaced (absolute -> log(t) is equidist.)
     */

    void run(double Tmax=25000, double deltaT=0.1, bool solve_implicit=true,
             size_t wint=500, size_t write_log=0) {
        fstream out(fn_concentration.c_str(), std::ios_base::out | std::ios_base::app);
        out.precision(25);

        size_t t0 = time(NULL);

        // Init solver and start
        vector_type x(initial_con), last_con(initial_con);
        double last_write(initial_t);

        // Function that is called by stepper and writes concentration (<vec>) to file.
        auto write_state = [this, &out, t0, wint, deltaT, Tmax, &last_con, &last_write]( const vector_type &vec , const double t ) {
            // calculate mean square distance from last write
            double last_msd=0;
            for(size_t i=0; i<vec.size(); ++i)
                 last_msd += (vec[i]-last_con[i])*(vec[i]-last_con[i]);

            if(t != last_write)
                last_msd /= (vec.size()*(t-last_write));
            last_write=t;
            last_con=vec;

            // write time and msd ...
            out << t << "," << last_msd;
            // and concentrations
            for(size_t l=0; l<vec.size(); ++l)
                out << "," << vec(l);

            out << std::endl;
        };

        size_t step_no = 0;  // Number of steps done by integrator (for diagnostics)
        // Calculate vector of time points at which "write_state" will be called.
        std::vector<double> times( wint );
        for( size_t i=0 ; i<wint ; ++i )
            if(write_log == 0)
                // linearly spaced output times
                times[i] = initial_t + double(i+1)/(wint)*(Tmax-initial_t);
            else if(write_log == 1 || initial_t == 0)
                // logarithmically spaced output times (relative to t0)
                times[i] = initial_t + (exp(double(i+1))-1)/(exp(double(wint))-1)*(Tmax-initial_t);
            else
                // if t0!=0 and one want's logarithmically spaced times (relative to t=0)!
                times[i] = initial_t*pow(exp(1/double(wint)*log(Tmax/initial_t)),double(i+1));

        // Call integrator - returns number of steps. Either implicit or explicit stepper is used.
        step_no = solve_implicit ?
                      integrate_times( make_controlled< rosenbrock4< double > >( 1.0e-6 , 1.0e-6 ) ,
                                       make_pair( stiff_system(*this) , stiff_system_jacobi(*this) ) ,
                                       x , times, deltaT, write_state) :
                      integrate_times( make_controlled< runge_kutta_dopri5< vector_type > >( 1.0e-6 , 1.0e-6 ) ,
                                       stiff_system(*this) , x, times, deltaT, write_state);

        // print time + steps needed and close file
        size_t t1 = time(NULL);
        std::cout << "Run took " << t1-t0 << " seconds and " << step_no << " steps!" << std::endl;
        out.close();
   }

   /*
    * The fastint integration (calculating forward and backward rates directly
    * through reaction constants) needs the networks to be exanded what this
    * method ensures (multiple occurance of same educt product in same reaction
    * has to be explicit - "A + A" instead of "2 A").
    */

    void initialize_fastint() {
        for(size_t i=0; i<re.size(); ++i) {
            for(size_t j=0; j<re[i].get_no_educt_s(); ++j) {
                if(re[i].get_educt_mul(j) > 1.1) {
                    re[i].add_educt_s(re[i].get_educt_id(j), re[i].get_educt_mul(j)-1);
                    re[i].set_educt_mul(j,1.0);
                }
            }

            for(size_t j=0; j<re[i].get_no_product_s(); ++j) {
                if(re[i].get_product_mul(j) > 1.1) {
                    re[i].add_product_s(re[i].get_product_id(j), re[i].get_product_mul(j)-1);
                    re[i].set_product_mul(j,1.0);
                }
            }
        }
    }

    /*
     * Method for calling the simple / fast integrator that uses explicit stepping and
     * the reaction constants defined in the network object. System (f(x)) is defined
     * in "fast_system".
     * For description of method's parameters see method "run" above.
     */

    void run_fastint(double Tmax=25000, double deltaT=0.1, size_t wint=500, size_t write_log=0) {
        fstream out(fn_concentration.c_str(), std::ios_base::out | std::ios_base::app);
        out.precision(25);

        size_t t0 = time(NULL);

        vector_type x(initial_con), last_con(initial_con);
        double last_write(initial_t);

        // Function that is called by stepper and writes concentration (<vec>) to file.
        auto write_state = [this, &out, t0, wint, deltaT, Tmax, &last_con, &last_write]( const vector_type &vec , const double t ) {
            // calculate mean square distance from last write
            double last_msd=0;
            for(size_t i=0; i<vec.size(); ++i)
                 last_msd += (vec[i]-last_con[i])*(vec[i]-last_con[i]);
            if(t != last_write)
                last_msd /= (vec.size()*(t-last_write));
            last_write=t;
            last_con=vec;

            // write time + msd...
            out << t << "," << last_msd;
            // ...and all concentrations
            for(size_t l=0; l<vec.size(); ++l)
                out << "," << vec(l);

            // If msd (change of concentration) per time and per species is less
            // than 1e-20 but not zero this integrator stops early because the
            // assumption is that the network is converged. (Might not be true in
            // the strict sense for all imaginable networks that can be simulated!)
            if(last_msd > 1e-44 && last_msd < 1e-20) {
                std::cout << "Reached msd < 1e-20 condition - exiting early!" << std::endl;
                size_t t1 = time(NULL);
                std::cout << "Run took " << t1-t0 << " seconds!" << std::endl;
                exit(0);
            }

            out << std::endl;
        };

        // Calculate vector of time points at which "write_state" will be called.
        // For details see same code in method "run" above.
        std::vector<double> times( wint );
        for( size_t i=0 ; i<wint ; ++i )
            if(write_log == 0)
                times[i] = initial_t + double(i+1)/(wint)*(Tmax-initial_t);
            else if(write_log == 1 || initial_t == 0)
                times[i] = initial_t + (exp(double(i+1))-1)/(exp(double(wint))-1)*(Tmax-initial_t);
            else
                times[i] = initial_t*pow(exp(1/double(wint)*log(Tmax/initial_t)),double(i+1));

        // Call integrator - returns number of steps
        size_t step_no = integrate_times( make_dense_output< runge_kutta_dopri5< vector_type > >( 1.0e-6 , 1.0e-6 ) ,
                         fast_system(*this) , x, times, deltaT, write_state);

        // print time + steps needed and close file
        size_t t1 = time(NULL);
        std::cout << "Run took " << t1-t0 << " seconds and " << step_no << " steps!" << std::endl;
        out.close();
   }
};

/*
 * Main-Function. Manages commandline parameter and prints help screen.
 * All the real work is done through using the reaction_network_system
 * class (above).
 */

int main(int argc, const char *argv []){
    srand(time(0));
    std::cout << "odeint_rnet version " << odeint_rnet_version_string << " (commit:"
              << GIT_VERSION << ")" << std::endl;

    cl_para cl(argc, argv);  // create class to query command line parameters

    // Prints the right hand side of the ODE (f(x)) for diagnostic purposes.

    if(cl.have_param("print_rhs")) {
        std::string fn_network, fn_concentration;

        if(cl.have_param("net"))
            fn_network=cl.get_param("net");
        else {
            cout << "You have to give the name of reaction network by 'net'!" << endl;
            return 1;
        }

        if(cl.have_param("con"))
            fn_concentration=cl.get_param("con");
        else {
            cout << "You have to give the name of the concentration file by 'con'!" << endl;
            return 1;
        }

        // load network and concentration file and call function to print right hand side
        reaction_network_system rns = reaction_network_system(fn_network, fn_concentration);
        rns.print_rhs();
    }

    // Simulates / integrates a reaction network while holding the boundary point
    // species (constant=true) constant. Effective reaction rates are calculated
    // directly from chemical energies. Stiff solver is available with option
    // "solve_implicit". Very slow for larger networks (>50 species).

    if(cl.have_param("simulate")) {
        std::string fn_network, fn_concentration;

        if(cl.have_param("net"))
            fn_network=cl.get_param("net");
        else {
            cout << "You have to give the name of reaction network by 'net'!" << endl;
            return 1;
        }

        if(cl.have_param("con"))
            fn_concentration=cl.get_param("con");
        else {
            cout << "You have to give the name of the concentration file by 'con'!" << endl;
            return 1;
        }

        // Initial time-step for integrator
        double deltaT = cl.have_param("deltaT") ? cl.get_param_d("deltaT") : 0.1;
        // ODE is integrated up to <Tmax>
        double Tmax = cl.have_param("Tmax") ? cl.get_param_d("Tmax") : 25000;
        // Number of times the output is written between initial time and <Tmax>
        double wint = cl.have_param("wint") ? cl.get_param_d("wint") : 500;
        // Use implicit solver?
        bool solve_implicit=cl.have_param("solve_implicit");
        // Is output given logarithmically spaced? (Or linearly?)
        size_t write_log=0;
        if(cl.have_param("write_log"))
            write_log=1;
        if(cl.have_param("write_log_abs"))
            write_log=2;

        std::cout << "Parameters are deltaT=" << deltaT << "  and Tmax=" << Tmax
                  << "   wint=" << wint << std::endl;

        if(write_log == 1)
            cout << "Period of output will be equidistant on logscale (from t=t_0)." << endl;

        if(write_log == 2)
            cout << "Period of output will be equidistant on logscale (from t=0)." << endl;

        if(solve_implicit)
            cout << "Stiff solver is used!" << endl;

        // Load reaction network and concentration file
        reaction_network_system rns = reaction_network_system(fn_network, fn_concentration);

        // Simulate ODE (and write results to file)
        // The method gives diagnostic feedback to the user through console output
        rns.run(Tmax, deltaT, solve_implicit, wint, write_log);
    }

    // Simulates / integrates a reaction network while holding the boundary point
    // species (constant=true) constant. Official parameter to call this is
    // "fastint", the option to call using "simsim" is maintained for backward
    // compatibility. This integrator is faster as the one above, especially for
    // larger (>50 species) networks. It uses the reaction constants present in
    // the network and no stiff solver.

    if(cl.have_param("fastint") || cl.have_param("simsim")) {
        std::string fn_network, fn_concentration;

        if(cl.have_param("net"))
            fn_network=cl.get_param("net");
        else {
            cout << "You have to give the name of reaction network by 'net'!" << endl;
            return 1;
        }

        if(cl.have_param("con"))
            fn_concentration=cl.get_param("con");
        else {
            cout << "You have to give the name of the concentration file by 'con'!" << endl;
            return 1;
        }

        // Initial time-step for integrator
        double deltaT = cl.have_param("deltaT") ? cl.get_param_d("deltaT") : 0.1;
        // ODE is integrated up to <Tmax>
        double Tmax = cl.have_param("Tmax") ? cl.get_param_d("Tmax") : 25000;
        // Number of times the output is written between initial time and <Tmax>
        double wint = cl.have_param("wint") ? cl.get_param_d("wint") : 500;
        // Is output given logarithmically spaced? (Or linearly?)
        size_t write_log=0;
        if(cl.have_param("write_log"))
            write_log=1;
        if(cl.have_param("write_log_abs"))
            write_log=2;

        std::cout << "Parameters are deltaT=" << deltaT << "  and Tmax=" << Tmax
                  << "   wint=" << wint << std::endl;

        if(write_log == 1)
            cout << "Period of output will be equidistant on logscale (from t=t_0)." << endl;

        if(write_log == 2)
            cout << "Period of output will be equidistant on logscale (from t=0)." << endl;

        // Load reaction network and concentration file
        reaction_network_system rns = reaction_network_system(fn_network, fn_concentration);

        // Simulate ODE (and write results to file)
        // The method gives diagnostic feedback to the user through console output
        rns.initialize_fastint();
        rns.run_fastint(Tmax, deltaT, wint, write_log);
    }

    // User interface / help dialogue
    if(cl.have_param("help") || cl.have_param("info")) {
        cout << "          odeint_rnet" << endl;
        cout << "          ===========" << endl;
        cout << " call with parameter 'info' or 'help' for showing this screen" << endl;
        cout << "current version is " << odeint_rnet_version_string << " (commit:" << GIT_VERSION << ")" << endl;
        cout << endl;
        cout << "-'print_rhs': load reaction network 'net', (last) concentration in 'con'! and" << endl;
        cout << "prints  right hand side of ode for diagnostic purposes." << endl;
	cout << endl;
        cout << "-'simulate': load reaction network 'net' and simulate file 'con'!. Parameters are:" << endl;
        cout << "   x'deltaT': Integration interval (relevant for integrator, not for output!)" << endl;
        cout << "   x'Tmax': Time up to which the system is simulated" << endl;
        cout << "   x'wint': number of output times between Tstart and Tmax" << endl;
        cout << "   x'write_log': write output logarithmically spaced (from t=t_0)" << endl;
        cout << "   x'write_log_abs': write output logarithmically spaced from (t=0)" << endl;
        cout << "   x'solve_implicit': use stiff solver (might be slower for big nets!)" << endl;
        cout << endl;
        cout << "-'fastint': load reaction network 'net' and simulate file 'con'!. Other than above" << endl;
        cout << "   this integrator does not offer a stiff solver. It also does not use energies to" << endl;
        cout << "   calculate effective rates directly but uses reaction constants. Parameters are:" << endl;
        cout << "   x'deltaT': Integration interval (relevant for integrator, not for output!)" << endl;
        cout << "   x'Tmax': Time up to which the system is simulated" << endl;
        cout << "   x'wint': number of output times between Tstart and Tmax" << endl;
        cout << "   x'write_log': write output logarithmically spaced (from t=t_0)" << endl;
        cout << "   x'write_log_abs': write output logarithmically spaced from (t=0)" << endl;
	cout << endl;
    }
}

ext: cpp


text:
#define BOOST_TEST_MAIN
#define BOOST_TEST_DYN_LINK

#include <memory>
#include <string>
#include <tuple>
#include <boost/test/unit_test.hpp>

#include "jml/arch/futex.h"
#include "jml/arch/timers.h"
#include "jml/utils/testing/watchdog.h"
#include "soa/service/message_loop.h"
#include "soa/service/http_client.h"
#include "soa/utils/print_utils.h"

#include "test_http_services.h"

using namespace std;
using namespace Datacratic;

/* helpers functions used in tests */
namespace {

typedef tuple<HttpClientError, int, string> ClientResponse;

#define CALL_MEMBER_FN(object, pointer)  (object.*(pointer))

/* sync request helpers */
template<typename Func>
ClientResponse
doRequest(MessageLoop & loop,
          const string & baseUrl, const string & resource,
          Func func,
          const RestParams & queryParams, const RestParams & headers,
          int timeout = -1)
{
    ClientResponse response;

    auto client = make_shared<HttpClient>(baseUrl, 4);
    loop.addSource("client", client);
    client->waitConnectionState(AsyncEventSource::CONNECTED);

    int done(false);
    auto onResponse = [&] (const HttpRequest & rq,
                           HttpClientError error,
                           int status,
                           string && headers,
                           string && body) {
        int & code = get<1>(response);
        code = status;
        string & body_ = get<2>(response);
        body_ = move(body);
        HttpClientError & errorCode = get<0>(response);
        errorCode = error;
        done = true;
        ML::futex_wake(done);
    };
    auto cbs = make_shared<HttpClientSimpleCallbacks>(onResponse);

    CALL_MEMBER_FN(*client, func)(resource, cbs, queryParams, headers,
                                  timeout);

    while (!done) {
        int oldDone = done;
        ML::futex_wait(done, oldDone);
    }

    loop.removeSource(client.get());
    client->waitConnectionState(AsyncEventSource::DISCONNECTED);

    return response;
}

ClientResponse
doGetRequest(MessageLoop & loop,
             const string & baseUrl, const string & resource,
             const RestParams & queryParams = RestParams(),
             const RestParams & headers = RestParams(),
             int timeout = -1)
{
    return doRequest(loop, baseUrl, resource, &HttpClient::get,
                     queryParams, headers, timeout);
}

ClientResponse
doDeleteRequest(MessageLoop & loop,
                const string & baseUrl, const string & resource,
                const RestParams & queryParams = RestParams(),
                const RestParams & headers = RestParams(),
                int timeout = -1)
{
    return doRequest(loop, baseUrl, resource, &HttpClient::del,
                     queryParams, headers, timeout);
}

ClientResponse
doUploadRequest(MessageLoop & loop,
                bool isPut,
                const string & baseUrl, const string & resource,
                const string & body, const string & type)
{
    ClientResponse response;

    auto client = make_shared<HttpClient>(baseUrl, 4);
    loop.addSource("client", client);
    client->waitConnectionState(AsyncEventSource::CONNECTED);

    int done(false);
    auto onResponse = [&] (const HttpRequest & rq,
                           HttpClientError error,
                           int status,
                           string && headers,
                           string && body) {
        int & code = get<1>(response);
        code = status;
        string & body_ = get<2>(response);
        body_ = move(body);
        HttpClientError & errorCode = get<0>(response);
        errorCode = error;
        done = true;
        ML::futex_wake(done);
    };

    auto cbs = make_shared<HttpClientSimpleCallbacks>(onResponse);
    HttpRequest::Content content(body, type);
    if (isPut) {
        client->put(resource, cbs, content);
    }
    else {
        client->post(resource, cbs, content);
    }

    while (!done) {
        int oldDone = done;
        ML::futex_wait(done, oldDone);
    }

    loop.removeSource(client.get());
    client->waitConnectionState(AsyncEventSource::DISCONNECTED);

    return response;
}

}

#if 1
BOOST_AUTO_TEST_CASE( test_http_client_get )
{
    cerr << "client_get\n";
    ML::Watchdog watchdog(10);
    auto proxies = make_shared<ServiceProxies>();
    HttpGetService service(proxies);

    service.addResponse("GET", "/coucou", 200, "coucou");
    service.start();

    MessageLoop loop;
    loop.start();

    service.waitListening();

#if 0
    /* request to bad ip
       Note: if the ip resolution timeout is very high on the router, the
       Watchdog timeout might trigger first */
    {
        ::fprintf(stderr, "request to bad ip\n");
        string baseUrl("http://123.234.12.23");
        auto resp = doGetRequest(loop, baseUrl, "/");
        BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::CouldNotConnect);
        BOOST_CHECK_EQUAL(get<1>(resp), 0);
    }
#endif

#if 0
    /* request to bad hostname
       Note: will fail when the name service returns a "default" value for all
       non resolved hosts */
    {
        ::fprintf(stderr, "request to bad hostname\n");
        string baseUrl("http://somewhere.lost");
        auto resp = doGetRequest(loop, baseUrl, "/");
        BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::HostNotFound);
        BOOST_CHECK_EQUAL(get<1>(resp), 0);
    }
#endif

    /* request with timeout */
    {
        ::fprintf(stderr, "request with timeout\n");
        string baseUrl("http://127.0.0.1:" + to_string(service.port()));
        auto resp = doGetRequest(loop, baseUrl, "/timeout", {}, {}, 1);
        BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::Timeout);
        BOOST_CHECK_EQUAL(get<1>(resp), 0);
    }

    /* request connection close */
    {
        ::fprintf(stderr, "testing behaviour with connection: close\n");
        string baseUrl("http://127.0.0.1:" + to_string(service.port()));
        auto resp = doGetRequest(loop, baseUrl, "/connection-close");
        BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::None);
        BOOST_CHECK_EQUAL(get<1>(resp), 204);
    }

    /* request to /nothing -> 404 */
    {
        ::fprintf(stderr, "request with 404\n");
        string baseUrl("http://127.0.0.1:"
                       + to_string(service.port()));
        auto resp = doGetRequest(loop, baseUrl, "/nothing");
        BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::None);
        BOOST_CHECK_EQUAL(get<1>(resp), 404);
    }

    /* request to /coucou -> 200 + "coucou" */
    {
        ::fprintf(stderr, "request with 200\n");
        string baseUrl("http://127.0.0.1:"
                       + to_string(service.port()));
        auto resp = doGetRequest(loop, baseUrl, "/coucou");
        BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::None);
        BOOST_CHECK_EQUAL(get<1>(resp), 200);
        BOOST_CHECK_EQUAL(get<2>(resp), "coucou");
    }

    /* headers and cookies */
    {
        string baseUrl("http://127.0.0.1:" + to_string(service.port()));
        auto resp = doGetRequest(loop, baseUrl, "/headers", {},
                                 {{"someheader", "somevalue"}});
        Json::Value expBody;
        expBody["accept"] = "*/*";
        expBody["host"] = baseUrl.substr(7);
        expBody["someheader"] = "somevalue";
        Json::Value jsonBody = Json::parse(get<2>(resp));
        BOOST_CHECK_EQUAL(jsonBody, expBody);
    }

    /* query-params */
    {
        string baseUrl("http://127.0.0.1:" + to_string(service.port()));
        auto resp = doGetRequest(loop, baseUrl, "/query-params",
                                 {{"value", "hello"}});
        string body = get<2>(resp);
        BOOST_CHECK_EQUAL(body, "?value=hello");
    }

    service.shutdown();
}
#endif

#if 1
BOOST_AUTO_TEST_CASE( test_http_client_post )
{
    cerr << "client_post\n";
    ML::Watchdog watchdog(10);
    auto proxies = make_shared<ServiceProxies>();
    HttpUploadService service(proxies);
    service.start();

    MessageLoop loop;
    loop.start();

    /* request to /coucou -> 200 + "coucou" */
    {
        string baseUrl("http://127.0.0.1:"
                       + to_string(service.port()));
        auto resp = doUploadRequest(loop, false, baseUrl, "/post-test",
                                    "post body", "application/x-nothing");
        BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::None);
        BOOST_CHECK_EQUAL(get<1>(resp), 200);
        Json::Value jsonBody = Json::parse(get<2>(resp));
        BOOST_CHECK_EQUAL(jsonBody["verb"], "POST");
        BOOST_CHECK_EQUAL(jsonBody["payload"], "post body");
        BOOST_CHECK_EQUAL(jsonBody["type"], "application/x-nothing");
    }

    service.shutdown();
}
#endif

#if 1
BOOST_AUTO_TEST_CASE( test_http_client_put )
{
    cerr << "client_put\n";
    ML::Watchdog watchdog(10);
    auto proxies = make_shared<ServiceProxies>();
    HttpUploadService service(proxies);
    service.start();

    MessageLoop loop;
    loop.start();

    string baseUrl("http://127.0.0.1:"
                   + to_string(service.port()));
    string bigBody;
    for (int i = 0; i < 65535; i++) {
        bigBody += "this is one big body,";
    }
    auto resp = doUploadRequest(loop, true, baseUrl, "/put-test",
                                bigBody, "application/x-nothing");
    BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::None);
    BOOST_CHECK_EQUAL(get<1>(resp), 200);
    Json::Value jsonBody = Json::parse(get<2>(resp));
    BOOST_CHECK_EQUAL(jsonBody["verb"], "PUT");
    BOOST_CHECK_EQUAL(jsonBody["payload"], bigBody);
    BOOST_CHECK_EQUAL(jsonBody["type"], "application/x-nothing");

    service.shutdown();
}
#endif

#if 1
BOOST_AUTO_TEST_CASE( http_test_client_delete )
{
    cerr << "client_delete" << endl;
    ML::Watchdog watchdog(10);

    auto proxies = make_shared<ServiceProxies>();
    HttpGetService service(proxies);

    service.addResponse("DELETE", "/deleteMe", 200, "Deleted");
    service.start();

    MessageLoop loop;
    loop.start();

    string baseUrl("http://127.0.0.1:" + to_string(service.port()));
    auto resp = doDeleteRequest(loop, baseUrl, "/deleteMe", {}, {}, 1);

    BOOST_CHECK_EQUAL(get<0>(resp), HttpClientError::None);
    BOOST_CHECK_EQUAL(get<1>(resp), 200);

    service.shutdown();
}
#endif

#if 1
BOOST_AUTO_TEST_CASE( test_http_client_put_multi )
{
    cerr << "client_put_multi\n";
    auto proxies = make_shared<ServiceProxies>();
    HttpUploadService service(proxies);
    service.start();

    string baseUrl("http://127.0.0.1:"
                   + to_string(service.port()));

    MessageLoop loop;
    loop.start();

    auto client = make_shared<HttpClient>(baseUrl);
    loop.addSource("client", client);
    client->waitConnectionState(AsyncEventSource::CONNECTED);

    size_t maxRequests(500);
    int done(0);

    auto makeBody = [&] (size_t i) {
        int multiplier = (i < maxRequests / 2) ? -2 : 2;
        size_t bodySize = 2000 + multiplier * i;
        string body = ML::format("%.4x", bodySize);
        size_t rndSize = bodySize - body.size();
        body += randomString(rndSize);

        return body;
    };

    for (size_t i = 0; i < maxRequests; i++) {
        auto sendBody = makeBody(i);
        auto onResponse = [&, sendBody] (const HttpRequest & rq,
                                         HttpClientError error,
                                         int status,
                                         string && headers,
                                         string && body) {
            BOOST_CHECK_EQUAL(error, HttpClientError::None);
            BOOST_CHECK_EQUAL(status, 200);
            Json::Value jsonBody = Json::parse(body);
            BOOST_CHECK_EQUAL(jsonBody["verb"], "PUT");
            BOOST_CHECK_EQUAL(jsonBody["payload"], sendBody);
            BOOST_CHECK_EQUAL(jsonBody["type"], "text/plain");
            done++;
            if (done == maxRequests) {
                ML::futex_wake(done);
            }
        };

        auto cbs = make_shared<HttpClientSimpleCallbacks>(onResponse);
        HttpRequest::Content content(sendBody, "text/plain");
        while (!client->put("/", cbs, content)) {
            ML::sleep(0.2);
        }
    };

    while (done < maxRequests) {
        int oldDone = done;
        ML::futex_wait(done, oldDone);
    }

    loop.shutdown();
    service.shutdown();
}
#endif

#if 1
/* Ensures that all requests are correctly performed under load, including
   when "Connection: close" is encountered once in a while.
   Not a performance test. */
BOOST_AUTO_TEST_CASE( test_http_client_stress_test )
{
    cerr << "stress_test\n";
    // const int mask = 0x3ff; /* mask to use for displaying counts */
    // ML::Watchdog watchdog(300);
    auto proxies = make_shared<ServiceProxies>();
    auto doStressTest = [&] (int numParallel) {
        ::fprintf(stderr, "stress test with %d parallel connections\n",
                  numParallel);

        HttpGetService service(proxies);
        service.start();
        service.waitListening();

        string baseUrl("http://127.0.0.1:"
                       + to_string(service.port()));

        MessageLoop loop;
        loop.start();

        auto client = make_shared<HttpClient>(baseUrl, numParallel);
        loop.addSource("client", client);
        client->waitConnectionState(AsyncEventSource::CONNECTED);

        int maxReqs(30000), numReqs(0), missedReqs(0);
        int numResponses(0);

        auto onDone = [&] (const HttpRequest & rq,
                           HttpClientError errorCode, int status,
                           string && headers, string && body) {
            numResponses++;

            BOOST_CHECK_EQUAL(errorCode, HttpClientError::None);
            BOOST_CHECK_EQUAL(status, 200);

            if (numResponses == numReqs) {
                ML::futex_wake(numResponses);
            }
        };

        while (numReqs < maxReqs) {
            const char * url = "/counter";
            auto cbs = make_shared<HttpClientSimpleCallbacks>(onDone);
            if (client->get(url, cbs)) {
                numReqs++;
                // if ((numReqs & mask) == 0 || numReqs == maxReqs) {
                //     ::fprintf(stderr, "performed %d requests\n", numReqs);
                // }
            }
            else {
                missedReqs++;
            }
        }

        ::fprintf(stderr, "all requests performed, awaiting responses...\n");
        while (numResponses < maxReqs) {
            int old(numResponses);
            ML::futex_wait(numResponses, old);
        }
        ::fprintf(stderr, "performed %d requests; missed: %d\n",
                  maxReqs, missedReqs);

        loop.shutdown();
        service.shutdown();
    };

    doStressTest(1);
    doStressTest(8);
    doStressTest(128);
}
#endif

#if 1
/* Ensure that the move constructor and assignment operator behave
   reasonably well. */
BOOST_AUTO_TEST_CASE( test_http_client_move_constructor )
{
    cerr << "move_constructor\n";
    ML::Watchdog watchdog(30);
    auto proxies = make_shared<ServiceProxies>();

    HttpGetService service(proxies);
    service.addResponse("GET", "/", 200, "coucou");
    service.start();
    service.waitListening();

    MessageLoop loop;
    loop.start();

    string baseUrl("http://127.0.0.1:"
                   + to_string(service.port()));

    auto doGet = [&] (HttpClient & getClient) {
        loop.addSource("client", getClient);
        getClient.waitConnectionState(AsyncEventSource::CONNECTED);

        int done(false);

        auto onDone = [&] (const HttpRequest & rq,
                           HttpClientError errorCode, int status,
                           string && headers, string && body) {
            done = true;
            ML::futex_wake(done);
        };
        auto cbs = make_shared<HttpClientSimpleCallbacks>(onDone);

        getClient.get("/", cbs);
        while (!done) {
            int old = done;
            ML::futex_wait(done, old);
        }

        loop.removeSource(&getClient);
        getClient.waitConnectionState(AsyncEventSource::DISCONNECTED);
    };

    /* move constructor */
    cerr << "testing move constructor\n";
    auto makeClient = [&] () {
        return HttpClient(baseUrl, 1);
    };
    HttpClient client1(move(makeClient()));
    doGet(client1);

    /* move assignment operator */
    cerr << "testing move assignment op.\n";
    HttpClient client2("http://nowhere", 1);
    client2 = move(client1);
    doGet(client2);

    service.shutdown();
}
#endif

#if 1
/* Ensure that an infinite number of requests can be queued when queue size is
 * 0, even from within callbacks. */
BOOST_AUTO_TEST_CASE( test_http_client_unlimited_queue )
{
    static const int maxLevel(4);

    ML::Watchdog watchdog(30);
    auto proxies = make_shared<ServiceProxies>();

    HttpGetService service(proxies);
    service.addResponse("GET", "/", 200, "coucou");
    service.start();
    service.waitListening();

    MessageLoop loop;
    loop.start();

    string baseUrl("http://127.0.0.1:"
                   + to_string(service.port()));

    auto client = make_shared<HttpClient>(baseUrl, 4, 0);
    loop.addSource("client", client);
    client->waitConnectionState(AsyncEventSource::CONNECTED);

    atomic<int> pending(0);
    int done(0);

    function<void(int)> doGet = [&] (int level) {
        pending++;
        auto onDone = [&,level] (const HttpRequest & rq,
                                 HttpClientError errorCode, int status,
                                 string && headers, string && body) {
            if (level < maxLevel) {
                for (int i = 0; i < 10; i++) {
                    doGet(level+1);
                }
            }
            pending--;
            done++;
        };
        auto cbs = make_shared<HttpClientSimpleCallbacks>(onDone);
        client->get("/", cbs);
    };

    doGet(0);

    while (pending > 0) {
        ML::sleep(1);
        cerr << "requests done: " + to_string(done) + "\n";
    }

    loop.shutdown();
    service.shutdown();
}
#endif

#if 1
/* Test connection restoration after a timeout occurs. */
BOOST_AUTO_TEST_CASE( test_http_client_connection_timeout )
{
    ML::Watchdog watchdog(30);
    auto proxies = make_shared<ServiceProxies>();

    HttpGetService service(proxies);
    service.addResponse("GET", "/", 200, "coucou");
    service.start();
    service.waitListening();

    MessageLoop loop;
    loop.start();

    string baseUrl("http://127.0.0.1:" + to_string(service.port()));

    auto client = make_shared<HttpClient>(baseUrl, 1);
    client->enableDebug(true);
    loop.addSource("client", client);
    client->waitConnectionState(AsyncEventSource::CONNECTED);

    int done(0);
    auto onDone = [&] (const HttpRequest & rq,
                       HttpClientError errorCode, int status,
                       string && headers, string && body) {
        done++;
        ML::futex_wake(done);
    };
    auto cbs = make_shared<HttpClientSimpleCallbacks>(onDone);
    client->get("/timeout", cbs, {}, {}, 1);
    client->get("/", cbs, {}, {}, 1);

    while (done < 2) {
        ML::futex_wait(done, done);
    }

    loop.shutdown();
    service.shutdown();
}
#endif

#if 1
/* Test connection restoration after the server closes the connection, under
 * various circumstances. */
BOOST_AUTO_TEST_CASE( test_http_client_connection_closed )
{
    ML::Watchdog watchdog(30);
    auto proxies = make_shared<ServiceProxies>();

    HttpGetService service(proxies);
    service.addResponse("GET", "/", 200, "coucou");
    service.start();
    service.waitListening();

    MessageLoop loop;
    loop.start();

    string baseUrl("http://127.0.0.1:" + to_string(service.port()));

    /* response sent, "Connection: close" header */
    {
        cerr << "* connection-close\n";
        auto client = make_shared<HttpClient>(baseUrl, 1);
        loop.addSource("client", client);
        client->waitConnectionState(AsyncEventSource::CONNECTED);

        int done(0);
        auto onDone = [&] (const HttpRequest & rq,
                           HttpClientError errorCode, int status,
                           string && headers, string && body) {
            done++;
            ML::futex_wake(done);
        };
        auto cbs = make_shared<HttpClientSimpleCallbacks>(onDone);
        client->get("/connection-close", cbs);
        client->get("/", cbs);

        while (done < 2) {
            ML::futex_wait(done, done);
        }

        loop.removeSourceSync(client.get());
    }

    /* Response sent, no "Connection: close" header. Performed multiple times
     * to ensure that extranenous close events are ignored appropriately when
     * the response was received and the actual connection released.
     * Unstable test due to the very specific condition during which the error
     * occurs: no request must be queued when the request is done processing
     * but one must be present when the closing of the connection is being
     * handled. */
    cerr << "* no connection-close\n";
    for (int i = 0; i < 10; i++)
    {
        auto client = make_shared<HttpClient>(baseUrl, 1);
        loop.addSource("client", client);
        client->waitConnectionState(AsyncEventSource::CONNECTED);

        atomic<int> done(0);
        auto onDone = [&] (const HttpRequest & rq,
                           HttpClientError errorCode, int status,
                           string && headers, string && body) {
            done++;
            ML::futex_wake(done);
        };
        auto cbs = make_shared<HttpClientSimpleCallbacks>(onDone);
        for (size_t i = 0; i < 3; i++) {
            client->post("/quiet-connection-close", cbs, string("no data"));
            while (done < i) {
                ML::futex_wait(done, done);
            }
        }

        ML::sleep(0.5);
        loop.removeSourceSync(client.get());
    }

    /* response not sent */
    {
        cerr << "* no response at all\n";
        auto client = make_shared<HttpClient>(baseUrl, 1);
        loop.addSource("client", client);
        client->waitConnectionState(AsyncEventSource::CONNECTED);

        int done(0);
        auto onDone = [&] (const HttpRequest & rq,
                           HttpClientError errorCode, int status,
                           string && headers, string && body) {
            done++;
            ML::futex_wake(done);
        };
        auto cbs = make_shared<HttpClientSimpleCallbacks>(onDone);
        client->get("/abrupt-connection-close", cbs);
        client->get("/", cbs);

        while (done < 2) {
            ML::futex_wait(done, done);
        }

        loop.removeSourceSync(client.get());
    }

    loop.shutdown();
    service.shutdown();
}
#endif

ext: cc

max_stars_repo_name: datacratic/soa

max_issues_repo_name: datacratic/soa

max_forks_repo_name: datacratic/soa


text:
/*
    Copyright (c) 2014, Philipp Krähenbühl
    All rights reserved.

    Redistribution and use in source and binary forms, with or without
    modification, are permitted provided that the following conditions are met:
        * Redistributions of source code must retain the above copyright
        notice, this list of conditions and the following disclaimer.
        * Redistributions in binary form must reproduce the above copyright
        notice, this list of conditions and the following disclaimer in the
        documentation and/or other materials provided with the distribution.
        * Neither the name of the Stanford University nor the
        names of its contributors may be used to endorse or promote products
        derived from this software without specific prior written permission.

    THIS SOFTWARE IS PROVIDED BY Philipp Krähenbühl ''AS IS'' AND ANY
    EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
    WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
    DISCLAIMED. IN NO EVENT SHALL Philipp Krähenbühl BE LIABLE FOR ANY
    DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
    (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
	 LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
	 ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
	 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
    SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/
#include "imgproc/color.h"
#include "imgproc/filter.h"
#include "imgproc/gradient.h"
#include "imgproc/morph.h"
#include "imgproc/nms.h"
#include "imgproc/resample.h"
#include "util/util.h"
#include "imgproc.h"
#include "gop.h"
#include "util.h"
#include <boost/python/def_visitor.hpp>

/****************  color.cpp  ****************/
typedef void (*ConvFunction)(Image &, const Image &);
static ConvFunction cfun[] = {rgb2luv,srgb2luv,rgb2lab,srgb2lab,rgb2hsv};
template<int T>
static Image convert( const Image & image ) {
	if( image.C() != 3 )
		throw std::invalid_argument( "3-channel image required" );
	Image r( image.W(), image.H(), 3 );
	cfun[T]( r, image );
	return r;
}

/****************  filter.cpp  ****************/
static RMatrixXf percentileFilter_m( const RMatrixXf & image, int rad, float p ) {
	RMatrixXf r( image.rows(), image.cols() );
	percentileFilter( r.data(), image.data(), image.cols(), image.rows(), 1, rad, p );
	return r;
}
#define defFilter( name )\
static RMatrixXf name##_m( const RMatrixXf & image, int rad ) {\
	RMatrixXf r( image.rows(), image.cols() );\
	name( r.data(), image.data(), image.cols(), image.rows(), 1, rad );\
	return r;\
}
defFilter( boxFilter )
defFilter( tentFilter )
defFilter( gaussianFilter )
static RMatrixXf exactGaussianFilter_m( const RMatrixXf & image, int rad, int R ) {
	RMatrixXf r( image.rows(), image.cols() );
	exactGaussianFilter( r.data(), image.data(), image.cols(), image.rows(), 1, rad, R );
	return r;
}

/****************  gradient.cpp  ****************/
static tuple gradient2( const Image & image ) {
	Image gx, gy;
	gradient( gx, gy, image );
	return make_tuple(gx, gy);
}
static tuple gradientMagAndOri2( const Image & image, int norm_rad, float norm_const ) {
	RMatrixXf gm, go;
	gradientMagAndOri( gm, go, image, norm_rad, norm_const );
	return make_tuple(gm, go);
}
static Image gradientHist2( const RMatrixXf & gm, const RMatrixXf & go, int nori, int nbins=1) {
	Image r;
	gradientHist( r, gm, go, nori, nbins );
	return r;
}

/****************  morph.cpp  ****************/
static RMatrixXb thin( const RMatrixXb & image ) {
	RMatrixXb r = image;
	thinningGuoHall( r );
	return r;
}

/****************  nms.cpp  ****************/
static RMatrixXf suppressBnd2( const RMatrixXf &im, int b ) {
	RMatrixXf r = im;
	suppressBnd( r, b );
	return r;
}
/****************  other  ****************/
static np::ndarray extractPatch( const np::ndarray &a, const RMatrixXi & l, const VectorXi & s ) {
	// Check the input
	if( !(a.get_flags() & np::ndarray::C_CONTIGUOUS) )
		throw std::invalid_argument( "Contiguous array required!" );
	if( l.cols() != a.get_nd() )
		throw std::invalid_argument("location dimension does not match matrix");
	if( l.cols() != s.size() )
		throw std::invalid_argument("patch size and locations have different dimension");
	for( int i=0; i<l.rows(); i++ )
		for( int j=0; j<s.size(); j++ )
			if( l(i,j)<0 || l(i,j)+s[j] > a.shape(j) )
				throw std::invalid_argument("location out of range! Dimension "+std::to_string(j)+" : "+std::to_string(l(i,j))+" + "+std::to_string(s[j])+"  ( 0 .. "+std::to_string(a.shape(j))+" )");

	const int itemsize = a.get_dtype().get_itemsize();

	// Create the output array
	std::vector<Py_intptr_t> out_size;
	out_size.push_back( l.rows() );
	for( int i=0; i<s.size(); i++ )
		if( s[i]>1 )
			out_size.push_back( s[i] );
	np::ndarray res = np::empty( out_size.size(), out_size.data(), a.get_dtype() );

	// Prepare the patches
	const int ls = s[s.size()-1]; // Radius of the last dimension
	std::vector<int> patch_offsets(1,0);
	for( int d=0; d<s.size()-1; d++ ) {
		std::vector<int> opo = patch_offsets;
		patch_offsets.clear();
		for( int i=0; i<s[d]; i++ )
			for( int o: opo )
				patch_offsets.push_back( a.shape(d+1)*(o+i) );
	}
	const char * pa = a.get_data();
	char * pres = res.get_data();
	for( int i=0,k=0; i<l.rows(); i++ ) {
		int o=0;
		for( int d=0; d<s.size(); d++ )
			o = a.shape(d)*o + l(i,d);
		// Extract the patches
		for( int po: patch_offsets ) {
			memcpy( pres+k*itemsize, pa+(po+o)*itemsize, ls*itemsize );
			k += ls;
		}
	}
	return res;
}
static np::ndarray extractBoxes( const Image &im, const RMatrixXi & b, int W, int H ) {
	// Create the output array
	const int C = im.C();
	const int im_W = im.W();

	np::ndarray res = np::empty( make_tuple( b.rows(), H, W, C ), np::dtype::get_builtin<float>() );
	const float * pim = (const float *)im.data();
	float * pres = (float *)res.get_data();
#pragma omp parallel for
	for( int i=0; i<b.rows(); i++ ) {
		// Build the lerp lookup table
		const float Y0=b(i,0), X0=b(i,1);
		const float dy=b(i,2)-Y0-1, dx=b(i,3)-X0-1;
		std::vector< int > x0( W ), y0( H ), x1( W ), y1( H );
		std::vector< float > wx( W ), wy( H );
		for( int j=0; j<H; j++ ) {
			float p = Y0 + j*dy/(H-1);
			y0[j] = floor(p); y1[j] = ceil(p);
			wy[j] = y1[j]-p;
		}
		for( int j=0; j<W; j++ ) {
			float p = X0 + j*dx/(W-1);
			x0[j] = floor(p); x1[j] = ceil(p);
			wx[j] = x1[j]-p;
		}

		// Run the lerp
		float * pr = pres + i*H*W*C;
		for( int j=0; j<H; j++ )
			for( int i=0; i<W; i++ )
				for( int k=0; k<C; k++ )
					pr[ (j*W+i)*C+k ] =    wx[i] *   wy[j] *pim[ ( y0[j]*im_W+x0[i] )*C+k ]
					                   +   wx[i] *(1-wy[j])*pim[ ( y1[j]*im_W+x0[i] )*C+k ]
					                   +(1-wx[i])*   wy[j] *pim[ ( y0[j]*im_W+x1[i] )*C+k ]
					                   +(1-wx[i])*(1-wy[j])*pim[ ( y1[j]*im_W+x1[i] )*C+k ];
	}
	return res;
}
template<typename T> struct ImageTypeStr {
	static const std::string s;
};
template<> const std::string ImageTypeStr<uint8_t>::s = "u1";
template<> const std::string ImageTypeStr<float>::s = "f4";

template<typename IM>
struct Image_indexing_suite: def_visitor<Image_indexing_suite<IM> > {
	typedef typename IM::value_type value_type;

	struct Image_pickle_suite : pickle_suite
	{
		static tuple getinitargs(const IM& im) {
			return make_tuple( im.W(), im.H(), im.C() );
		}

		static object getstate(const IM& im) {
			const int N = im.W()*im.H()*im.C()*sizeof(value_type);
			return object( handle<>( PyBytes_FromStringAndSize( (const char*)im.data(), N ) ) );
		}

		static void setstate(IM& im, const object & state) {
			if(!PyBytes_Check(state.ptr()))
				throw std::invalid_argument("Failed to unpickle, unexpected type!");
			const int N = im.W()*im.H()*im.C()*sizeof(value_type);
			if( PyBytes_Size(state.ptr()) != N )
				throw std::invalid_argument("Failed to unpickle, unexpected size!");
			memcpy( im.data(), PyBytes_AS_STRING(state.ptr()), N );
		}
	};

	template <class classT>	void visit(classT& c) const {
		c
		.def("__init__",make_constructor(&Image_indexing_suite::init1))
		.def("__init__",make_constructor(&Image_indexing_suite::init2))
		.def("__init__",make_constructor(&Image_indexing_suite::init3))
		.add_property("W",&IM::W)
		.add_property("H",&IM::H)
		.add_property("C",&IM::C)
		.def("tileC",&IM::tileC)
		.def_pickle(Image_pickle_suite())
		.add_property("__array_interface__", &Image_indexing_suite::array_interface);
	}
	static IM * init1() {
		return new IM();
	}
	static IM * init2( int W, int H, int C ) {
		return new IM( W, H, C );
	}
	static IM * init3( const np::ndarray & d ) {
		checkArray( d, value_type, 2, 3, true );

		IM* r = new IM(d.shape(1),d.shape(0),d.get_nd()>2?d.shape(2):1);
		memcpy( r->data(), d.get_data(), r->W()*r->H()*r->C()*sizeof(value_type) );
		return r;
	}
	static dict array_interface( IM & im ) {
		dict r;
		r["shape"] = make_tuple( im.H(), im.W(), im.C() );
		r["typestr"] = ImageTypeStr<value_type>::s;
		r["data"] = make_tuple((size_t)im.data(),1);
		r["version"] = 3;
		return r;
	}
};
template<typename I1, typename I2>
I1 convertImage( const I2 & im ) {
	return im;
}
void defineImgProc() {
	ADD_MODULE(imgproc);
	// Color
	def("rgb2luv",convert<0>);
	def("srgb2luv",convert<1>);
	def("rgb2lab",convert<2>);
	def("srgb2lab",convert<3>);
	def("rgb2hsv",convert<4>);
	// Filters
	def("boxFilter",boxFilter_m);
	def("tentFilter",tentFilter_m);
	def("gaussianFilter",gaussianFilter_m);
	def("exactGaussianFilter",exactGaussianFilter_m);
	def("percentileFilter",percentileFilter_m);
	// Gradient
	def("gradient",gradient2);
	def("gradientMag",gradientMag);
	def("gradientMagAndOri",gradientMagAndOri2);
	def("gradientHist",gradientHist2);
	// Morphology
	def("thin",thin);
	// NMS
	def("nms",nms);
	def("suppressBnd",suppressBnd2);
	// Upsampling
	def("upsample",upsample);
	def("upsampleLinear",upsampleLinear);
	def("downsample", (RMatrixXf(*)( const RMatrixXf &, int, int ))downsample);
	def("downsample", (Image(*)( const Image &, int, int ))downsample);
	def("padIm",padIm);
	// Other
	def("extractPatch",extractPatch);
	def("extractBoxes",extractBoxes);

	def("imread",imread);
	def("imwrite",imwrite);

	// Image and Image8u
	class_<Image,std::shared_ptr<Image> >("Image")
	.def(Image_indexing_suite<Image>())
	.def("toImage8u",convertImage<Image,Image8u>);
	class_<Image8u,std::shared_ptr<Image8u> >("Image8u")
	.def(Image_indexing_suite<Image8u>())
	.def("toImage",convertImage<Image,Image8u>);
}

max_stars_repo_name: PingCheng-Wei/SD-MaskRCNN

max_issues_repo_name: PingCheng-Wei/SD-MaskRCNN

max_forks_repo_name: PingCheng-Wei/SD-MaskRCNN


text:
/**
 * @file process/parsers/SearchParser.cpp
 * @author Ian Yang
 * @date Created <2010-06-11 17:12:27>
 */
#include "SearchParser.h"

#include <common/BundleSchemaHelpers.h>
#include <common/Keys.h>

#include <boost/algorithm/string/case_conv.hpp>

namespace sf1r {

using driver::Keys;

/**
 * @class SearchParser
 *
 * The @b search field is an object. It specifies a full text search in
 * documents. The search object has following fields:
 *
 * - @b keywords* (@c String): Keywords used in the search.
 * - @b USERID (@c String): USERID is used to get the user info, which assists
 *   to get personalized search result. This item is not required to set. if it
 *   is not set, the result is the normal ranking search result.
 * - @b in (@c Array): Which properties the search is performed in. Every item
 *   can be an Object or an String. If it is an Object, the \b property field is
 *   used as the property name, otherwise, the String itself is used as property
 *   name. If this is not specified, all indexed properties are used. Valid
 *   properties can be check though schema/get (see SchemaController::get() ).
 *   All index names (the @b name field in every index) can be used.
 * - @b count (@c Array): Which properties the COUNT() is performed in. Every item
 *   can be an Object or an String. If it is an Object, the \b property field is
 *   used as the property name, otherwise, the String itself is used as property
 *   name. Valid properties should have numeric type.
 * - @b group_label (@c Array): Get documents in the specified group labels. Each label consists of a property name and value.@n
 *   You could specify multiple labels, for example, there are 4 labels named as "A1", "A2", "B1", "B2".@n
 *   If the property name of label "A1" and "A2" is "A", while the property name of label "B1" and "B2" is "B",@n
 *   then the documents returned would belong to the labels with the condition of ("A1" or "A2") and ("B1" or "B2").
 *   - @b property* (@c String): the property name.
 *   - @b value* (@c Array): the label value.@n
 *   For the property type of string, it's an array of the path from root to leaf node.
 *   Each is a @c String of node value. Only those documents belonging to the leaf node
 *   would be returned.@n
 *   For the property type of int or float, you could either specify a property value
 *   or a range in @b value[0]. If it is a range, you could not specify the same @b property
 *   in request["group"].@n
 *   Regarding the range format, you could specify it in the form of "101-200",
 *   meaning all values between 101 and 200, including both boundary values,
 *   or in the form of "-200", meaning all values not larger than 200,
 *   or in the form of "101-", meaning all values not less than 101.@n
 *   For the property type of datetime, in @b value[0], you could specify a string of year, month or day.@n
 *   For example, "2012" means the year of 2012, "2012-07" means the month of 2012 July, "2012-07-06" means the day of 2012 July 6th.
 *   - @b auto_select_limit (@c Uint = 1): when @b value is an empty array, the top labels would be selected automatically.@n
 *   You could specify this parameter as how many top labels need to be selected automatically.@n
 *   The selected top labels are returned as @b top_group_label in search response.@n
 *   This feature only supports the property type of string.
 * - @b attr_label (@c Array): Get documents in the specified attribute labels. Each label consists of an attribute name and value.@n
 *   You could specify multiple labels, for example, there are 4 labels named as "A1", "A2", "B1", "B2".@n
 *   If the attribute name of label "A1" and "A2" is "A", while the attribute name of label "B1" and "B2" is "B",@n
 *   then the documents returned would belong to the labels with the condition of ("A1" or "A2") and ("B1" or "B2").
 *   - @b attr_name* (@c String): the attribute name
 *   - @b attr_value* (@c String): the attribute value
 * - @b ranking_model (@c String): How to rank the search result. Result with
 *   high ranking score is considered has high relevance. Now SF1 supports
 *   following ranking models.
 *   - @b plm
 *   - @b bm25
 *   - @b kl.@n
 *   If this is omitted, A default ranking model specified in configuration file
 *   is used (In Collection node attribute).
 * - @b searching_mode (@c Object): Searching mode options
 *   We can choose the searching mode to combine our
 *   query parsed term for supporting different search in different applications.
 *   following searching modes are available
 *   - @b and space among terms are explained as boolean AND
 *   - @b or space among terms are explained as boolean OR
 *   - @b wand wand means weak and, and \b threshold is used to indicate
 *   how close it is to boolean query: the closer \b threshold to 1, the closer it
 *   is to boolean AND
 *   - @b knn: simhash based query. The query is hashed into a fingerprint,
 *   while the search results are those documents whose fingerprints are within
 *   a threshold hamming distance
 *   - @b suffix  using suffix array based succint index to perform queries.
 *   \b lucky means top K number during this retrieval, using lucky means there
 *   exists approximation since they are not the real top K results returned due
 *   to performance issue.   \b use_fuzzy is a switch to indicate whether bag-of-words
 *   based fuzzy queries are enabled for this query. If it is set as false, it means
 *   the search results only contain the longest suffix of query.
 *   If this is omitted, @b and searching mode is used as the default value.
 *   - @b zambezi search in the zambezi index (an in-memory inverted index).
 * - @b log_keywords (@c Bool = @c true): Whether the keywords should be
 *   logged.
 * - @b analyzer (@c Object): Keywords analyzer options
 *   - @b apply_la (@c Bool = @c true): TODO
 *   - @b use_synonym_extension (@c Bool = @c false): TODO
 *   - @b use_original_keyword (@c Bool = @c false): TODO
 * - @b is_random_rank (@c Bool = @c false): If true, the search results would
 *   be randomly ordered.@n
 *   In order to enable this feature, you need also configure a non-zero weight
 *   for <ProductRanking><Score type="random"> in collection config file.
 * - @b is_require_related (@c Bool = @c false): If true, the search results would
 *   be contain related queries.
 * - @b query_source (@c String): Where does the query come from, used to decide
 *   the categories to boost in product ranking.
 * - @b boost_group_label (@c Array): The group labels to boost product rankings.@n
 *   It's an array of group labels, each element is an array of the label path
 *   started from root node. As an example, two labels could be written as
 *   [["手机数码", "手机通讯", "手机"], ["电脑办公", "电脑整机"]].@n
 *   In order to enable this feature, in collection config file
 *   <ProductRanking><Score type="category">, you need to configure the property
 *   name of these group labels, and also a non-zero weight.
 *
 * Example
 * @code
 * {
 *   "keywords": "test",
 *   "USERID": "56",
 *   "session_id":"session_123"
 *   "in": [
 *     {"property": "title"},
 *     {"property": "content"}
 *   ],
 *   "count": [
 *     {"property": "visitation"},
 *   ],
 *   "ranking_model": "plm",
 *   "searching_mode": {
 *     "mode": "wand",
 *     "threshold": "0.35"
 *    },
 *   "log_keywords": true,
 *   "is_require_related": true,
 *   "analyzer": {
 *     "use_synonym_extension": true,
 *     "apply_la": true,
 *     "use_original_keyword": true
 *   }
 * }
 * @endcode
 */

bool SearchParser::parse(const Value& search)
{
    clearMessages();

    // keywords
    keywords_ = asString(search[Keys::keywords]);
    if (keywords_.empty())
    {
        error() = "Require keywords in search.";
        return false;
    }

    userID_ = asString(search[Keys::USERID]);
    sessionID_ = asString(search[Keys::session_id]);
    querySource_ = asString(search[Keys::query_source]);

    int labelCount = 0;
    if (search.hasKey(Keys::taxonomy_label))
        ++labelCount;
    if (search.hasKey(Keys::name_entity_item) ||
         search.hasKey(Keys::name_entity_type))
        ++labelCount;
    if (labelCount > 1)
    {
        error() = "At most one label type (taxonomy label, name entity label) can be specided.";
        return false;
    }

    if (!parseGroupLabel_(search) ||
        !parseAttrLabel_(search) ||
        !parseAdSearch_(search) ||
        !parseBoostGroupLabel_(search))
        return false;

    logKeywords_ = asBoolOr(search[Keys::log_keywords], true);
    isRandomRank_ = asBoolOr(search[Keys::is_random_rank], false);
    requireRelatedQueries_ = asBoolOr(search[Keys::is_require_related], false);

    // counter list properties
    const Value& countNode = search[Keys::count];
    if (! nullValue(countNode))
    {
        if (countNode.type() == Value::kArrayType)
        {
            countList_.resize(countNode.size());
            for (std::size_t i = 0; i < countList_.size(); ++i)
            {
                const Value& currentProperty = countNode(i);
                if (currentProperty.type() == Value::kObjectType)
                {
                    countList_[i] = asString(currentProperty[Keys::property]);
                }
                else
                {
                    countList_[i] = asString(currentProperty);
                }

                if (countList_[i].empty())
                {
                    error() = "Failed to parse properties in count.";
                    return false;
                }

                // validation
                PropertyConfig propertyConfig;
                propertyConfig.setName(countList_[i]);
                if (!getPropertyConfig(indexSchema_, propertyConfig))
                {
                    error() = "Unknown property in count/in: " +
                              propertyConfig.getName();
                    return false;
                }

                if (!propertyConfig.bIndex_ ||!propertyConfig.bFilter_ || !propertyConfig.isNumericType())
                {
                    error() = "Counted property should be numeric filter type: " +
                                propertyConfig.getName();
                    return false;
                }
            }
        }
        else
        {
            error() = "Count must be an array";
            return false;
        }
    }

    // La
    const Value& analyzer = search[Keys::analyzer];
    analyzerInfo_.applyLA_ = asBoolOr(analyzer[Keys::apply_la], true);
    analyzerInfo_.useOriginalKeyword_ = asBool(analyzer[Keys::use_original_keyword]);
    analyzerInfo_.synonymExtension_ = asBool(analyzer[Keys::use_synonym_extension]);

    // ranker
    rankingModel_ = RankingType::DefaultTextRanker;
    if (search.hasKey(Keys::ranking_model))
    {
        Value::StringType ranker = asString(search[Keys::ranking_model]);
        boost::to_lower(ranker);

        if (ranker == "plm")
        {
            rankingModel_ = RankingType::PLM;
        }
        else if (ranker == "bm25")
        {
            rankingModel_ = RankingType::BM25;
        }
        else if (ranker == "kl")
        {
            rankingModel_ = RankingType::KL;
        }
        else
        {
            warning() = "Unknown rankingModel. Default one is used.";
        }
    }

    searchingModeInfo_.clear();

    //search mode
    const Value& searching_mode = search[Keys::searching_mode];
    if ( searching_mode.hasKey(Keys::mode) )
    {
        Value::StringType mode = asString(searching_mode[Keys::mode]);
        boost::to_lower(mode);
        if (mode == "and")
        {
            searchingModeInfo_.mode_ = SearchingMode::DefaultSearchingMode;
            if (searching_mode.hasKey(Keys::query_prune) && !asBool(searching_mode[Keys::query_prune]))
            {
                searchingModeInfo_.useQueryPrune_ = false;
            }
        }
        else if (mode == "or")
        {
            searchingModeInfo_.mode_ = SearchingMode::OR;
        }
        else if (mode == "wand")
        {
            searchingModeInfo_.mode_ = SearchingMode::WAND;
        }
        else if (mode == "verbose")
        {
            searchingModeInfo_.mode_ = SearchingMode::VERBOSE;
        }
        else if (mode == "suffix")
        {
            searchingModeInfo_.mode_ = SearchingMode::SUFFIX_MATCH;
            if (searching_mode.hasKey(Keys::query_prune) && !asBool(searching_mode[Keys::query_prune]))
            {
                searchingModeInfo_.useQueryPrune_ = false;
            }

            if (searching_mode.hasKey(Keys::use_fuzzyThreshold) && asBool(searching_mode[Keys::use_fuzzyThreshold]))
            {
                searchingModeInfo_.useFuzzyThreshold_ = true;
                if (searching_mode.hasKey(Keys::fuzzy_threshold))
                {
                    float threshold = (float)asDouble(searching_mode[Keys::fuzzy_threshold]);
                    if (threshold >= 0.1F && threshold <= 1.0F)
                    {
                        searchingModeInfo_.fuzzyThreshold_ = threshold;
                    }
                    else
                        warning() = "fuzzy threshold is invalid, must between 0.1 and 1.0, now set as Default 0.5";
                }

                if (searching_mode.hasKey(Keys::tokens_threshold))
                {
                    float threshold = (float)asDouble(searching_mode[Keys::tokens_threshold]);
                    if (threshold >= 0.1F && threshold <= 1.0F)
                    {
                        searchingModeInfo_.tokensThreshold_ = threshold;
                    }
                    else
                        warning() = "tokensThreshold_ is invalid, must between 0.1 and 1.0, now set as Default 0.5";
                }
            }

            if (searching_mode.hasKey(Keys::use_pivilegeQuery) && asBool(searching_mode[Keys::use_pivilegeQuery]))
            {
                searchingModeInfo_.usePivilegeQuery_ = true;
                if (searching_mode.hasKey(Keys::privilege_Query))
                {
                    std::string privilegeQuery = asString(searching_mode[Keys::privilege_Query]);
                    searchingModeInfo_.privilegeQuery_ =  privilegeQuery;
                }

                if (searching_mode.hasKey(Keys::privilege_Weight))
                {
                    float privilege_Weight = (float)asDouble(searching_mode[Keys::privilege_Weight]);
                    if (privilege_Weight >= 0.1F && privilege_Weight <= 1.0F)
                    {
                        searchingModeInfo_.privilegeWeight_ = privilege_Weight;
                    }
                    else
                        warning() = "privilegeWeight is invalid, must between 0.1 and 1.0, now set as Default 0.1";
                }
            }
        }
        else if (mode == "zambezi")
        {
            searchingModeInfo_.mode_ = SearchingMode::ZAMBEZI;

            // add Zambezi Algorithm

            if (searching_mode.hasKey(Keys::algorithm) )
            {
                Value::StringType algorithm = asString(searching_mode[Keys::algorithm]);
                boost::to_lower(algorithm);
                if (algorithm == "svs")
                {
                    searchingModeInfo_.algorithm_ = izenelib::ir::Zambezi::SVS;
                }
                else if (algorithm == "wand")
                {
                    searchingModeInfo_.algorithm_ = izenelib::ir::Zambezi::WAND;
                }
                else if (algorithm == "bwand_and")
                {
                    searchingModeInfo_.algorithm_ = izenelib::ir::Zambezi::BWAND_AND;
                }
                else if (algorithm == "bwand_or")
                {
                    searchingModeInfo_.algorithm_ = izenelib::ir::Zambezi::BWAND_OR;
                }
                else if (algorithm == "mbwand")
                {
                    searchingModeInfo_.algorithm_ = izenelib::ir::Zambezi::MBWAND;
                }
                else
                {
                    error() = "The Zambezi search algorithm:" + algorithm + " is wrong";
                    return false;
                }
            }
        }
        else if (mode == "ad_index")
        {
            searchingModeInfo_.mode_ = SearchingMode::AD_INDEX;
        }
        else if (mode == "ad_recommend")
        {
            searchingModeInfo_.mode_ = SearchingMode::AD_RECOMMEND;
        }
        else if (mode == "sponsored_ad_search")
        {
            searchingModeInfo_.mode_ = SearchingMode::SPONSORED_AD_SEARCH;
        }
        else
        {
            warning() = "Unknown searchingMode. Default searching mode is used.";
        }

        if (indexSchema_.empty() && searchingModeInfo_.mode_ !=SearchingMode::ZAMBEZI)
        {
            error() = "IndexSchema is not in xml config, the search mode is wrong";
            return false;
        }

        if (searching_mode.hasKey(Keys::threshold))
        {
            float threshold = (float)asDouble(searching_mode[Keys::threshold]);
            if (threshold < 1)
            {
                searchingModeInfo_.threshold_ = threshold;
            }
            else
            {
                warning() = "Threshold is invalid. Warning: threshold must be smaller than one";
            }
        }

        if (searching_mode.hasKey(Keys::lucky))
        {
            uint32_t lucky = asUint(searching_mode[Keys::lucky]);
            if (lucky <= 100000)
            {
                searchingModeInfo_.lucky_ = lucky;
            }
            else
            {
                warning() = "Threshold is invalid. Warning: lucky must be not greater than than 1000";
            }
        }

        if (searching_mode.hasKey(Keys::original_query))
        {
            searchingModeInfo_.useOriginalQuery_ = asBool(searching_mode[Keys::original_query]);
        }

        if (searching_mode.hasKey(Keys::use_fuzzy))
        {
            searchingModeInfo_.usefuzzy_ = asBool(searching_mode[Keys::use_fuzzy]);
        }
        if (searching_mode.hasKey(Keys::filter_mode))
        {
            Value::StringType filter_mode = asString(searching_mode[Keys::filter_mode]);
            boost::to_lower(filter_mode);
            if(filter_mode == "or")
            {
                searchingModeInfo_.filtermode_ = SearchingMode::OR_Filter;
            }
            else if(filter_mode == "and")
            {
                searchingModeInfo_.filtermode_ = SearchingMode::AND_Filter;
            }
            else
            {
                warning() = "Unknown filter Mode. Default filtering mode is used.";
            }
        }

    }

    // properties
    const Value& propertiesNode = search[Keys::in];
    if (nullValue(propertiesNode))
    {
        if (searchingModeInfo_.mode_ != SearchingMode::ZAMBEZI)
            getDefaultSearchPropertyNames(indexSchema_, properties_);
        else
            getDefaultZambeziSearchPropertyNames(zambeziConfig_, properties_);

    }
    else if (propertiesNode.type() == Value::kArrayType)
    {
        properties_.resize(propertiesNode.size());
        for (std::size_t i = 0; i < properties_.size(); ++i)
        {
            const Value& currentProperty = propertiesNode(i);
            if (currentProperty.type() == Value::kObjectType)
            {
                properties_[i] = asString(currentProperty[Keys::property]);
            }
            else
            {
                properties_[i] = asString(currentProperty);
            }

            if (properties_[i].empty())
            {
                error() = "Failed to parse properties in search.";
                return false;
            }

            // validation
            PropertyConfig propertyConfig;
            propertyConfig.setName(properties_[i]);

            if (searchingModeInfo_.mode_ == SearchingMode::ZAMBEZI)
            {
                if (!zambeziConfig_.checkValidationIndexConfig(properties_[i]))
                {
                    error() = "Unknown property in search/in Using zambezi: " +
                      propertyConfig.getName();
                    return false;
                }
                continue;
            }
            else if (!getPropertyConfig(indexSchema_, propertyConfig) )
            {
                error() = "Unknown property in search/in: " +
                      propertyConfig.getName();
                return false;
            }

            if (!(propertyConfig.bIndex_||propertyConfig.bSuffixIndex_))
            {
                warning() = "Property is not indexed, ignore it: " +
                            propertyConfig.getName();
            }
        }
    }

    if (properties_.empty())
    {
        error() = "Require list of properties in which search is performed.";
        return false;
    }

    return true;
}

bool SearchParser::parseGroupLabel_(const Value& search)
{
    const Value& groupLabelArray = search[Keys::group_label];

    if (nullValue(groupLabelArray))
        return true;

    if (groupLabelArray.type() != Value::kArrayType)
    {
        error() = "Require an array for group labels.";
        return false;
    }

    for (std::size_t i = 0; i < groupLabelArray.size(); ++i)
    {
        const Value& groupLabelElem = groupLabelArray(i);
        std::string propName = asString(groupLabelElem[Keys::property]);

        if (propName.empty())
        {
            error() = "Must specify both property name and array of value path for group label";
            return false;
        }

        const Value& pathNode = groupLabelElem[Keys::value];
        faceted::GroupParam::GroupPath groupPath;
        if (pathNode.type() == Value::kArrayType)
        {
            for (std::size_t j = 0; j < pathNode.size(); ++j)
            {
                std::string nodeValue = asString(pathNode(j));
                if (nodeValue.empty())
                {
                    error() = "request[search][group_label][value] is empty for property " + propName;
                    return false;
                }
                groupPath.push_back(nodeValue);
            }
        }
        else if (pathNode.type() != Value::kNullType)
        {
            error() = "Require an array for request[search][group_label][value]";
            return false;
        }

        if (groupPath.empty())
        {
            int limit = asUintOr(groupLabelElem[Keys::auto_select_limit], 1);
            if (limit == 0)
            {
                error() = "Must specify a positive value for 'auto_select_limit'";
                return false;
            }
            groupLabelAutoSelectLimits_[propName] = limit;
        }
        else
        {
            groupLabels_[propName].push_back(groupPath);
        }
    }

    return true;
}

bool SearchParser::parseAttrLabel_(const Value& search)
{
    const Value& attrNode = search[Keys::attr_label];

    if (nullValue(attrNode))
        return true;

    if (attrNode.type() != Value::kArrayType)
    {
        error() = "Require an array for attr labels.";
        return false;
    }

    for (std::size_t i = 0; i < attrNode.size(); ++i)
    {
        const Value& attrPair = attrNode(i);
        std::string attrName = asString(attrPair[Keys::attr_name]);
        std::string attrValue = asString(attrPair[Keys::attr_value]);

        if (attrName.empty() || attrValue.empty())
        {
            error() = "Must specify both attribute name and value for attr label";
            return false;
        }

        attrLabels_[attrName].push_back(attrValue);
    }

    return true;
}

bool SearchParser::parseAdSearch_(const Value& search)
{
    const Value& adSearchNode = search[Keys::ad_search];

    if (nullValue(adSearchNode))
        return true;

    if (adSearchNode.type() != Value::kArrayType)
    {
        error() = "Require an array for attr labels.";
        return false;
    }

    for(std::size_t i = 0; i< adSearchNode.size(); i++)
    {
        const Value& valuePair = adSearchNode(i);
        std::string property = asString(valuePair[Keys::property]);
        std::string value = asString(valuePair[Keys::value]);

        if(property.empty() || value.empty())
        {
            error() = "Must specify both attribute name and value ";
            return false;
        }
        adSearch_.push_back(std::make_pair(property, value));
    }

    return true;
}

bool SearchParser::parseBoostGroupLabel_(const Value& search)
{
    const Value& boostLabelValue = search[Keys::boost_group_label];

    if (nullValue(boostLabelValue))
        return true;

    if (boostLabelValue.type() != Value::kArrayType)
    {
        error() = "Require an array for search parameter boost_group_label.";
        return false;
    }

    for (std::size_t i = 0; i < boostLabelValue.size(); ++i)
    {
        const Value& pathValue = boostLabelValue(i);

        if (pathValue.type() != Value::kArrayType)
        {
            error() = "Require an array for each element in search parameter boost_group_label.";
            return false;
        }

        faceted::GroupParam::GroupPath groupPath;
        for (std::size_t j = 0; j < pathValue.size(); ++j)
        {
            std::string nodeStr = asString(pathValue(j));
            groupPath.push_back(nodeStr);
        }

        boostGroupLabels_.push_back(groupPath);
    }

    return true;
}

} // namespace sf1r