Linear Controller

#include "../flowstar/Continuous.h"
#include<fstream>
#include<ctime>

using namespace std;
using namespace flowstar;

int main(int argc, char* argv[])
{
    // Declaration of the state variables.
    unsigned int numVars = 3;

    int x0_id = stateVars.declareVar("x0");
    int x1_id = stateVars.declareVar("x1");
    int u_id = stateVars.declareVar("u");

    int domainDim = numVars + 1;

    // Define the continuous dynamics.
    Expression_AST<Real> deriv_x0("40 - x1");  // theta_r = 0
    Expression_AST<Real> deriv_x1("-0.2*x1 + u");
    Expression_AST<Real> deriv_u("0");

    vector<Expression_AST<Real>> ode_rhs(numVars);
    ode_rhs[x0_id] = deriv_x0;
    ode_rhs[x1_id] = deriv_x1;
    ode_rhs[u_id] = deriv_u;

    Deterministic_Continuous_Dynamics dynamics(ode_rhs);

    // Specify the parameters for reachability computation.
    Computational_Setting setting;

    unsigned int order = 12;
    setting.setFixedStepsize(0.01, order);
    setting.setTime(0.1);
    setting.setCutoffThreshold(1e-10);
    setting.setQueueSize(1000);
    setting.setRemainderEstimation(vector<Interval>(numVars, Interval(-0.01, 0.01)));
    setting.prepare();

    float state_1 = std::stof(argv[1]);
    float state_2 = std::stof(argv[2]);

    Interval init_x0(state_1 - 1, state_1 + 1), init_x1(state_2 - 2, state_2 + 2), init_u(0);
    vector<Interval> X0 = {init_x0, init_x1, init_u};

    Flowpipe initial_set(X0);

    // Result of the reachability computation
    Result_of_Reachability result;

    string theta_1 = argv[3];
    string theta_2 = argv[4];
    int step_num = std::stoi(argv[5]);

    for (int iter = 0; iter < step_num; ++iter)
    {
        vector<Interval> box;
        initial_set.intEval(box, order, setting.tm_setting.cutoff_threshold);

        string strExpU = "x0*" + theta_1 + "+x1*" + theta_2;

        Expression_AST<Real> exp_u(strExpU);

        TaylorModel<Real> tm_u;
        exp_u.evaluate(tm_u, initial_set.tmvPre.tms, order, initial_set.domain, setting.tm_setting.cutoff_threshold, setting.g_setting);

        initial_set.tmvPre.tms[u_id] = tm_u;

        dynamics.reach(result, setting, initial_set, vector<Constraint>());

        if (result.status == COMPLETED_SAFE || result.status == COMPLETED_UNSAFE || result.status == COMPLETED_UNKNOWN)
        {
            initial_set = result.fp_end_of_time;
        }
        else
        {
            printf("Terminated due to too large overestimation.\n");
            break;
        }
    }

    ofstream result_output("./outputs/reach_result.txt");
    vector<Interval> end_box;
    result.fp_end_of_time.intEval(end_box, order, setting.tm_setting.cutoff_threshold);
    string reach_result = "Verification result: Unknown(" + std::to_string(step_num) + ")";

    if (end_box[0].inf() >= -0.4 && end_box[0].sup() <= -0.28 && end_box[1].inf() >= 0.05 && end_box[1].sup() <= 0.22)
    {
        reach_result = "Verification result: Yes(" + std::to_string(step_num) + ")";
    }
    else if (end_box[0].inf() >= 0.28 || end_box[0].sup() <= -0.4 || end_box[1].inf() >= 0.22 || end_box[1].sup() <= 0.05)
    {
        reach_result = "Verification result: No(" + std::to_string(step_num) + ")";
    }

    if (result_output.is_open())
    {
        result_output << reach_result << endl;
    }

    return 0;
}