// source code for mixed-trust optimal E
// 
// Copyright 2024 Carnegie Mellon University.
// 
// Licensed under a -style license, please see license.txt or contact permission@sei.cmu.edu for full terms.
// 
// This material is based upon work funded and supported by the Department of Defense under Contract No. FA8702-15-D-0002 with Carnegie Mellon University for the operation of the Software Engineering Institute, a federally funded research and development center.  
// 
// The view, opinions, and/or findings contained in this material are those of the author(s) and should not be construed as an official Government position, policy, or decision, unless designated by other documentation. NO WARRANTY. THIS CARNEGIE MELLON UNIVERSITY AND SOFTWARE ENGINEERING INSTITUTE MATERIAL IS FURNISHED ON AN "AS-IS" BASIS. CARNEGIE MELLON UNIVERSITY MAKES NO WARRANTIES OF ANY KIND, EITHER EXPRESSED OR IMPLIED, AS TO ANY MATTER INCLUDING, BUT NOT LIMITED TO, WARRANTY OF FITNESS FOR PURPOSE OR MERCHANTABILITY, EXCLUSIVITY, OR RESULTS OBTAINED FROM USE OF THE MATERIAL. CARNEGIE MELLON UNIVERSITY DOES NOT MAKE ANY WARRANTY OF ANY KIND WITH RESPECT TO FREEDOM FROM PATENT, TRADEMARK, OR COPYRIGHT INFRINGEMENT.
// 
// [DISTRIBUTION STATEMENT A] This material has been approved for public release and unlimited distribution.  Please see Copyright notice for non-US Government use and distribution.
//
// DM24-0049
//
// mtoe.c
//
// What does this program do?
//   It performs optimal configuration for mixed-trust tasks.
//   Specifically, if there exists an assignment of E so that
//   the schedulability equations are satisfied, then the algorithm finds
//   an assignment so that the schedulability equations are satisfied.
//
// How do I install it?
//   Install z3. It will be in /usr/bin/z3. So we don't need to set the PATH variable.
//   Compile the program as follows:
//     gcc -o mtoe mtoe.c -O3 -lm
//   Now, you have a file mtoe. The file "mtoe" means Mixed-Trust Optimal E.
//		    
// How do I actually run it?
//   Open a console. Then type:
//      ./mtoe -i system.txt -o system_results_from_analysis.txt
//   With this command, you tell the tool that it should read the file system.txt and perform configuration and
//   write the output of the analysis to the file results_from_analysis.txt.
//   You can run it in meta mode (a mode where the program generates tasksets and configures them) with
//     ./mtoe -m &> outputlog.txt
//
// What assumptions are made?
//   It is assumed that each mixed-trust task has a hypertask; that is, we do not allow kC_i=0.
//
// What is the format of the input file?
//   The input file is a text file with numbers and these numbers are not comma separated. 
//   An example of an input file is as follows:
// 1
// 2
// 1 2
//   20 20 10 8 3
// 2 1
//   25 25 12 7 2
//
// The 1st number is "1"  and it indicates the method to be used.
// The 2nd number is "2"  and it indicates the number of tasks.

// The 3rd number is "1"  and it indicates that the following information is for task 1.
// The 4th number is "2"  and it indicates the priority of task 1.
// The 5th number is "20" and it indicates the T parameter of task 1.
// The 6th number is "20" and it indicates the D parameter of task 1.
// The 7th number is "10" and it indicates the E parameter of task 1 (this parameter is ignored for now).
// The 8th number is "8"  and it indicates the C parameter of task 1.
// The 9th number is "3"  and it indicates the kC parameter of task 1.
// Analogous for the remaining lines.
// Note that priorities must be unique.
//   
//  You can learn about the input file format by reading the function "void readtasksetfromfile(char* fn)" in this source code file.
//
// This source code is structured as follows:
//   1. Data that we use to get bounds on E
//   2. Basic routines, general
//   3. Basic routines for taskset
//   4. Code for parsing result file from Z3
//   5. Functions to compute response times of hypertasks
//   6. Functions to compute bounds on response times of guesttasks
//   7. Functions to add lines to SMT file
//   8. Functions for getting variable names
//   9. Functions for adding variables
//  10. Functions for adding constraints
//  11. Functions for retrieving information from result file from SMT solving
//  12. Functions for bounds on E
//  13. Functions for solving PR
//  14. Function for optimal configuration for assigning only E
//  15. Functions for generating random taskset
//  16. Main
//  17. Test results

#include <stdlib.h>
#include <stdio.h>
#include <time.h> // we add this in order to call the function time( ..)
#include <fcntl.h>
#include <math.h>
#include <sys/stat.h> // this is needed for the constant S_IWRITE
#include <string.h>   //  this is needed for strcmp
#include <unistd.h>
#include <stdio.h>
#include <stdlib.h>

#define MAXNTASKS      10000
#define MAXNVARSUSED   1000000

#define UTIL_EPSILON   0.001

#define VERBOSE 0 

int metathing = 0;

int outputprioritiesasnumbers = 1;

char inputfilename[20000]="system.txt";
char outputfilename[20000]="system_results_from_analysis.txt";

struct task {
  int id;
  int prio;
  double T;
  double D;
  double E;
  double C;
  double kC;
};

int usemethod;
int ntasks;
struct task tasks[MAXNTASKS+1];  // there is no element in index 0.

int   counterforSMTsolver = 0;
FILE* smt_instance_f; 
char  smt_instance_fn[200];
FILE* smt_cmd_f; 
char  smt_cmd_fn[200];
FILE* smt_result_f; 
char smt_result_fn[200];

//   1. Data that we use to get bounds on E

double LBE[MAXNTASKS+1];  // there is no element in index 0.
double UBE[MAXNTASKS+1];  // there is no element in index 0.
double LBRg[MAXNTASKS+1];  // there is no element in index 0.

//   2. Basic routines, general

double min2double(double t1,double t2) {
  if (t1<=t2) { return t1; }
  else        { return t2; } 
}
double max2double(double t1,double t2) {
  if (t1>=t2) { return t1; }
  else        { return t2; } 
}
double relu_ceildiv(double t1,double t2) {
  return max2double(ceil(t1/t2),0);
}
double ceildiv(double t1,double t2) {
  return ceil(t1/t2);
}
void initialize_random_generator() {
  srand(137);
  srand48(137);
}

int myrandom() { return rand(); }
double myrand48() { return drand48(); }
int my_int_random(int lb,int ub)             { return lb + (myrandom() % (ub-lb+1) ); }
double my_double_random(double lo,double hi) { return lo + myrand48()*(hi-lo);        }
double mypow(double t1,double t2) { return exp(t2*log(t1)); }
void setuserexecutablemode(char* fn) {
  char s[2000];
  sprintf(s,"chmod u+x %s", fn);
  system(s);
}
double log2(double t) { return (log(t))/(log(2.0)); }

//   3. Basic routines for taskset

double getutil() {
  int j; double sum;
  sum = 0.0;
  for (j=1;j<=ntasks;j++) {
    sum = sum + (tasks[j].C+tasks[j].kC)/tasks[j].T;
  }
  return sum;
}
int priorities_are_unique() {
  int j; int jprime;
  for (j=1;j<=ntasks-1;j++) { 
    for (jprime=j+1;jprime<=ntasks;jprime++) { 
      if (tasks[j].prio==tasks[jprime].prio) {
        return 0;      
      }
    }
  }
  return 1;
}
char get_x_from_integer(int an_integer) {
  if      (an_integer==0) { return 'E'; }
  else if (an_integer==1) { return 'A'; }
  else                    { printf("Error in get_x_from_integer. Incorrect value of an_integer.\n"); exit(-1); }
}

char get_y_from_integer(int an_integer) {
  if      (an_integer==0) { return 'E'; }
  else if (an_integer==1) { return 'A'; }
  else                    { printf("Error in get_y_from_integer. Incorrect value of an_integer.\n"); exit(-1); }
}
void readtasksetfromfile_task(FILE* myfile,int i) {
  fscanf(myfile,"%d %d",&(tasks[i].id),&(tasks[i].prio) );
  fscanf(myfile,"%lf %lf %lf %lf %lf",&(tasks[i].T),&(tasks[i].D),&(tasks[i].E),&(tasks[i].C),&(tasks[i].kC));
}
void readtasksetfromfile(char* fn) {
  FILE* myfile; int i;
  myfile = fopen( fn, "rt");
  if (myfile==NULL) { printf("Error opening file in function readtasksetfromfile. fn=%s\n",fn); fflush(stdout); exit(0); }
  fscanf(myfile,"%d",&usemethod);
  fscanf(myfile,"%d",&ntasks);
  for (i=1;i<=ntasks;i++) { readtasksetfromfile_task(myfile,i); }
  fclose(myfile);
}
void writetasksettofile_task(FILE* myfile, int i) {
  fprintf(myfile,"%d %d\n",tasks[i].id, tasks[i].prio );
  fprintf(myfile,"  %lf %lf %lf %lf %lf\n",tasks[i].T,tasks[i].D,tasks[i].E,tasks[i].C,tasks[i].kC );
}
void writetasksettofile(char* fn) {
  FILE* myfile; int i;
  myfile = fopen(fn,"w");
  if (myfile==NULL) { printf("Error opening file in function writetasksettofile. fn=%s\n", fn); fflush(stdout); exit(0); }
  fprintf( myfile, "%d\n", usemethod);
  fprintf( myfile, "%d\n", ntasks);
  for (i=1;i<=ntasks;i++) { writetasksettofile_task(myfile,i); }
  fclose( myfile);
}
#if COMPILEMSVS
unsigned __int128 diff_us_timeval(struct timeval *x, struct timeval *y)
{
  unsigned __int128 amillion = 1000000;
  unsigned __int128 u = ((x->tv_sec * amillion) + x->tv_usec);
  unsigned __int128 v = ((y->tv_sec * amillion) + y->tv_usec);
  return u-v;
}
double diff_s_timeval(struct timeval *x, struct timeval *y)
{
  unsigned __int128 amillion = 1000000;
  unsigned __int128 t;
  unsigned __int128 tsec;
  unsigned __int128 tusec;
  double t_as_double_secondpart;
  double t_as_double_usecondpart;
  double t_as_double;
  t = diff_us_timeval( x, y);
  tsec = t / amillion;
  tusec =  t % amillion; 
  t_as_double_secondpart = tsec;
  t_as_double_usecondpart = tusec;
  t_as_double_usecondpart = t_as_double_usecondpart / amillion;
  t_as_double = t_as_double_secondpart + t_as_double_usecondpart; 
  return t_as_double;
}
#else
unsigned long long diff_ns_timespec(struct timespec *x, struct timespec *y) {
  unsigned long long abillion = 1000000000;
  unsigned long long u;
  u = x->tv_sec;
  u = u * abillion;
  u = u + x->tv_nsec;

  unsigned long long v;
  v = y->tv_sec;
  v = v * abillion;
  v = v + y->tv_nsec;
  return u-v;
}
double diff_s_timespec (struct timespec *x, struct timespec *y)
{
  unsigned long long abillion = 1000000000;
  unsigned long long t;
  unsigned long long tsec;
  unsigned long long tnsec;
  double t_as_double_secondpart;
  double t_as_double_usecondpart;
  double t_as_double;
  t = diff_ns_timespec( x, y);
  tsec = t / abillion;
  tnsec =  t % abillion; 
  t_as_double_secondpart = tsec;
  t_as_double_usecondpart = tnsec;
  t_as_double_usecondpart = t_as_double_usecondpart / abillion;
  t_as_double = t_as_double_secondpart + t_as_double_usecondpart; 
  return t_as_double;
}
#endif

#if COMPILEMSVS

#define CLOCK_REALTIME 0
LARGE_INTEGER getFILETIMEoffset() {
  SYSTEMTIME s;
  FILETIME f;
  LARGE_INTEGER t;
  s.wYear = 1970;
  s.wMonth = 1;
  s.wDay = 1;
  s.wHour = 0;
  s.wMinute = 0;
  s.wSecond = 0;
  s.wMilliseconds = 0;
  SystemTimeToFileTime(&s, &f);
  t.QuadPart = f.dwHighDateTime;
  t.QuadPart <<= 32;
  t.QuadPart |= f.dwLowDateTime;
  return (t);
}
int my_clock_gettime(int X, struct timeval *tv) {
  LARGE_INTEGER           t;
  FILETIME            f;
  double                  microseconds;
  static LARGE_INTEGER    offset;
  static double           frequencyToMicroseconds;
  static int              initialized = 0;
  static BOOL             usePerformanceCounter = 0;
  if (!initialized) {
    LARGE_INTEGER performanceFrequency;
    initialized = 1;
    usePerformanceCounter = QueryPerformanceFrequency(&performanceFrequency);
    if (usePerformanceCounter) {
      QueryPerformanceCounter(&offset);
      frequencyToMicroseconds = (double)performanceFrequency.QuadPart / 1000000.;
    } else {
      offset = getFILETIMEoffset();
      frequencyToMicroseconds = 10.;
    }
  }
  if (usePerformanceCounter) { QueryPerformanceCounter(&t); }
  else {
    GetSystemTimeAsFileTime(&f);
    t.QuadPart = f.dwHighDateTime;
    t.QuadPart <<= 32;
    t.QuadPart |= f.dwLowDateTime;
  }
  t.QuadPart -= offset.QuadPart;
  microseconds = (double)t.QuadPart / frequencyToMicroseconds;
  t.QuadPart = microseconds;
  tv->tv_sec = t.QuadPart / 1000000;
  tv->tv_usec = t.QuadPart % 1000000;
  return (0);
}
#endif

//   4. Code for parsing result file from z3

int is_the_smt_result_satisfiable() {
  char s[200];
  smt_result_f = fopen(smt_result_fn,"rt");
  if (smt_result_f==NULL) {
    printf("In is_the_smt_result_satisfiable. Error when opening file. fn=%s\n", smt_result_fn); fflush(stdout);
    exit(0);
  }
  fgets(s,200,smt_result_f);
  fclose(smt_result_f); 
  if (strncmp(s,"sat",3)==0) {
    return 1;
  } else {
    return 0;
  }
}
int char_is_blank(char c) {
  if ((c==' ') || (c=='\t') || (c=='\n') || (c=='\r') || (c=='\f') || (c=='\v')) {
    return 1;  
  } else {
    return 0;   
  }
}
int trim_space_before(char* p) {
  int l; int i;
  l = strlen(p);
  if (char_is_blank(p[0])) {
    for (i=0;i<l-1;i++) {
      p[i] = p[i+1];    
    }
    p[l-1] = 0;
    return 1;
  } else {
    return 0;
  }
}
int trim_space_after(char* p) {
  int l;
  l = strlen(p);
  if (char_is_blank(p[l-1])) {
    p[l-1] = 0;
    return 1;
  } else {
    return 0;
  }
}
int trim_leftparenthesis_before(char* p) {
  int l; int i;
  l = strlen(p);
  if (p[0]=='(') {
    for (i=0;i<l-1;i++) {
      p[i] = p[i+1];    
    }
    p[l-1] = 0;
    return 1;
  } else {
    return 0;
  }
}
int trim_rightparenthesis_after(char* p) {
  int l;
  l = strlen(p);
  if (p[l-1]==')') {
    p[l-1] = 0;    
    return 1;
  } else {
    return 0; 
  }
}
int trim_spaces_and_parentheses_aux(char* p) {
  int changed;
  changed = 0;
  changed = changed || trim_space_before(p);
  changed = changed || trim_space_after(p);
  changed = changed || trim_leftparenthesis_before(p);
  changed = changed || trim_rightparenthesis_after(p);
  return changed;
}
void trim_spaces_and_parentheses(char* p) {
  int changed;
  changed = trim_spaces_and_parentheses_aux(p);
  while (changed) {
    changed = trim_spaces_and_parentheses_aux(p);
  }
}

//   3. Interact with human user

void parsearguments( int argc, char *argv[], char *envp[]) {
  int i;
  i = 1;
  while (i<=argc-1) {
    if        (strcmp( argv[i], "-i")==0) {
      strcpy( inputfilename, argv[i+1] );
      i = i + 2;
    } else if (strcmp( argv[i], "-o")==0) {
      strcpy( outputfilename, argv[i+1] );
      i = i + 2;
    } else if (strcmp( argv[i], "-m")==0) {
      metathing = 1;
      i = i + 1;
    } else {
      i++;
    }
  }
}
void printE_to_fn(char* fn) {
  FILE* f; int i;
  f = fopen(fn,"w");
  if (f==NULL) { printf("In printE_to_fn. fn=%s\n",fn); fflush(stdout); exit(0); }
  for (i=1;i<=ntasks;i++) {
    fprintf(f,"E_%d = %lf\n",i,tasks[i].E);
  }
  fclose(f);
}

//   5. Functions to compute response times of hypertasks

double get_max_kC_of_lp(int i) {
  double ret_val; int j;
  if (ntasks<=0) { printf("Error in get_max_kC_of_lp.\n"); exit(-1); }
  ret_val = 0.0;
  for (j=1;j<=ntasks;j++) { 
    if (tasks[j].prio<tasks[i].prio) {
      if (ret_val<tasks[j].kC) { ret_val=tasks[j].kC; }
    }
  }
  return ret_val;
}
double compute_rhs_of_eq3(int i,double t) {
  double term1; double term2; double term3; int j;
  term1 = get_max_kC_of_lp(i);
  term2 = ceil(t/tasks[i].T)*tasks[i].kC;
  term3 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio>tasks[i].prio) {
      term3 = term3 + ceil(t/tasks[j].T)*tasks[j].kC;
    }
  }
  return term1+term2+term3;
}
double compute_tiK(int i) {
  double t; double newt;
  if (tasks[i].kC<=0.0) {  printf("Error in compute_tiK\n"); exit(-1); }
  t = 0;
  newt = tasks[i].kC; 
  while (newt>t) {
    t = newt;
    newt = compute_rhs_of_eq3(i,t);
  }
  return t;
}
double compute_rhs_of_eq4(int i,int q,double w) {
  double term1; double term2; double term3; int j;
  term1 = get_max_kC_of_lp(i);
  term2 = (q-1)*tasks[i].kC;
  term3 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio>tasks[i].prio) {
      term3 = term3 + (floor(w/tasks[j].T)+1)*tasks[j].kC;
    }
  }
  return term1+term2+term3;
}
double compute_wiqK(int i,int q) {
  double w; double neww;
  w = 0.0;
  neww = compute_rhs_of_eq4(i,q,w);
  // note that if neww==0.0 here, that does not indicate an error; this is correct behavior if ntasks=1 and q=1.
  while (neww>w) {
    w = neww;
    neww = compute_rhs_of_eq4(i,q,w);
  }
  return w;
}
double compute_wiqK_and_RiK_and_return_RiK(int i,double tK) {
  double RiK; int q; double current_wiqK; double current_RiqK; 
  if (ceil(tK/tasks[i].T)<=0) { printf("Error in compute_wiqK_and_RiK_and_return_RiK.\n"); exit(-1); }
  RiK = 0.0;
  for (q=1;q<=ceil(tK/tasks[i].T);q++) {
    current_wiqK = compute_wiqK(i,q);
    current_RiqK = current_wiqK + tasks[i].kC - (q-1)*tasks[i].T;
    if (RiK<current_RiqK) {
      RiK = current_RiqK;    
    }
  }
  return RiK;
}
double compute_hypertask_responsetime(int i) {
  double t;
  t = compute_tiK(i);
  return compute_wiqK_and_RiK_and_return_RiK(i,t); 
}
int hypertask_response_times_respect_deadlines(double* RK) {
  int i;
  for (i=1;i<=ntasks;i++) {
    if (RK[i]>tasks[i].D) { return 0; }
  }
  return 1;  
}

//   6. Functions to compute bounds on response times of guesttasks

double LBrbf_evaluation(int j,char y,double t,int b) {
  double H1; double H2; double term1; double term2;
  if        (y=='E') {
    H1 = relu_ceildiv(t-(tasks[j].T-LBE[j]),tasks[j].T);
    H2 = ceildiv(t,tasks[j].T);
    term1 = H1*tasks[j].C*b;
    term2 = H2*tasks[j].kC;
    return term1+term2;
  } else if (y=='A') {
    H1 = ceildiv(t,tasks[j].T);
    H2 = relu_ceildiv(t-UBE[j],tasks[j].T);
    term1 = H1*tasks[j].C*b;
    term2 = H2*tasks[j].kC;
    return term1+term2;
  } else {
    printf("Error in LBrbf_evaluation\n"); exit(-1);  
  }
}
double UBrbf_evaluation(int j,char y,double t,int b) {
  double H1; double H2; double term1; double term2;
  if        (y=='E') {
    H1 = relu_ceildiv(t-(tasks[j].T-UBE[j]),tasks[j].T);
    H2 = ceildiv(t,tasks[j].T);
    term1 = H1*tasks[j].C*b;
    term2 = H2*tasks[j].kC;
    return term1+term2;
  } else if (y=='A') {
    H1 = ceildiv(t,tasks[j].T);
    H2 = relu_ceildiv(t-LBE[j],tasks[j].T);
    term1 = H1*tasks[j].C*b;
    term2 = H2*tasks[j].kC;
    return term1+term2;
  } else {
    printf("Error in UBrbf_evaluation\n"); exit(-1);  
  }
}
double compute_LB_rhs_of_eq7(int i,char x,double t) {
  double term1; double term2; double term3; int j;
  term1 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio<tasks[i].prio) {
      term1 = term1 + LBrbf_evaluation(j,'E',t,0);
    }
  }
  term2 = LBrbf_evaluation(i,x,t,1);
  term3 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio>tasks[i].prio) {
      term3 = term3 + max2double(LBrbf_evaluation(j,'E',t,1),LBrbf_evaluation(j,'A',t,1));
    }
  }
  return term1+term2+term3;
}
double compute_UB_rhs_of_eq7(int i,char x,double t) {
  double term1; double term2; double term3; int j;
  term1 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio<tasks[i].prio) {
      term1 = term1 + UBrbf_evaluation(j,'E',t,0);
    }
  }
  term2 = UBrbf_evaluation(i,x,t,1);
  term3 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio>tasks[i].prio) {
      term3 = term3 + max2double(UBrbf_evaluation(j,'E',t,1),UBrbf_evaluation(j,'A',t,1));
    }
  }
  return term1+term2+term3;
}
double compute_tLB(int i,char x) {
  double t; double newt;
  if (min2double(tasks[i].kC,tasks[i].C)<=0.0) {  printf("Error in compute_tLB\n"); exit(-1); }
  t = 0;
  newt = min2double(tasks[i].kC,tasks[i].C); 
  while (newt>t) {
    t = newt;
    newt = compute_LB_rhs_of_eq7(i,x,t);
  }
  return t;
}
double compute_tUB(int i,char x) {
  double t; double newt;
  if (min2double(tasks[i].kC,tasks[i].C)<=0.0) {  printf("Error in compute_tUB\n"); exit(-1); }
  t = 0;
  newt = min2double(tasks[i].kC,tasks[i].C); 
  while (newt>t) {
    t = newt;
    newt = compute_UB_rhs_of_eq7(i,x,t);
  }
  return t;
}
double compute_LB_rhs_of_eq8(int i,int q,char x,double w) {
  double term1; double term2; double term3; int j;
  term1 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio<tasks[i].prio) {
      term1 = term1 + LBrbf_evaluation(j,'E',w,0);
    }
  }
  if (x=='E') {
    term2 = q*tasks[i].C+q*tasks[i].kC;
  } else {
    term2 = q*tasks[i].C+(q-1)*tasks[i].kC;  
  }
  term3 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio>tasks[i].prio) {
      term3 = term3 + max2double(LBrbf_evaluation(j,'E',w,1),LBrbf_evaluation(j,'A',w,1));
    }
  }
  return term1+term2+term3;
}
double compute_UB_rhs_of_eq8(int i,int q,char x,double w) {
  double term1; double term2; double term3; int j;
  term1 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio<tasks[i].prio) {
      term1 = term1 + UBrbf_evaluation(j,'E',w,0);
    }
  }
  if (x=='E') {
    term2 = q*tasks[i].C+q*tasks[i].kC;
  } else {
    term2 = q*tasks[i].C+(q-1)*tasks[i].kC;  
  }
  term3 = 0.0;
  for (j=1;j<=ntasks;j++) {
    if (tasks[j].prio>tasks[i].prio) {
      term3 = term3 + max2double(UBrbf_evaluation(j,'E',w,1),UBrbf_evaluation(j,'A',w,1));
    }
  }
  return term1+term2+term3;
}
double compute_wjobLB(int i,int q,char x) {
  double w; double neww;
  if (tasks[i].C<=0.0) { printf("Error in compute_wjobLB\n"); exit(-1); }
  w = 0;
  if (x=='E') {
    neww = q*tasks[i].C+q*tasks[i].kC;
  } else {
    neww = q*tasks[i].C+(q-1)*tasks[i].kC;  
  }
  while (neww>w) {
    w = neww;
    neww = compute_LB_rhs_of_eq8(i,q,x,w);
  }
  return w;
}
double compute_wjobUB(int i,int q,char x) {
  double w; double neww;
  if (tasks[i].C<=0.0) { printf("Error in compute_wjobUB\n"); exit(-1); }
  w = 0;
  if (x=='E') {
    neww = q*tasks[i].C+q*tasks[i].kC;
  } else {
    neww = q*tasks[i].C+(q-1)*tasks[i].kC;  
  }
  while (neww>w) {
    w = neww;
    neww = compute_UB_rhs_of_eq8(i,q,x,w);
  }
  return w;
}
double compute_RgjobLB(int i,int q,char x) {
  double w; double RgjobLB;
  w = compute_wjobLB(i,q,x);
  if        (x=='E') {
    RgjobLB = w-((q-1)*tasks[i].T+(tasks[i].T-LBE[i])); 
  } else if (x=='A') {
    RgjobLB = w-(q-1)*tasks[i].T; 
  } else {
     printf("Error in compute_RgjobLB.\n"); exit(-1);  
  }
  return RgjobLB;
}
double compute_RgjobUB(int i,int q,char x) {
  double w; double RgjobUB;
  w = compute_wjobUB(i,q,x);
  if        (x=='E') {
    RgjobUB = w-((q-1)*tasks[i].T+(tasks[i].T-UBE[i])); 
  } else if (x=='A') {
    RgjobUB = w-(q-1)*tasks[i].T; 
  } else {
     printf("Error in compute_RgjobUB.\n"); exit(-1);  
  }
  return RgjobUB;
}
double compute_RgtaskxLB(int i,char x) {
  double RgtaskxLB; double t; int q; double currentjob_r;
  RgtaskxLB = 0.0;
  t = compute_tLB(i,x);
  for (q=1;q<=ceil(t/tasks[i].T);q++) {
    currentjob_r = compute_RgjobLB(i,q,x);
    if (RgtaskxLB<currentjob_r) { RgtaskxLB = currentjob_r; }
  }
  return RgtaskxLB;
}
double compute_RgtaskxUB(int i,char x) {
  double RgtaskxUB; double t; int q; double currentjob_r;
  RgtaskxUB = 0.0;
  t = compute_tUB(i,x);
  for (q=1;q<=ceil(t/tasks[i].T);q++) {
    currentjob_r = compute_RgjobUB(i,q,x);
    if (RgtaskxUB<currentjob_r) { RgtaskxUB = currentjob_r; }
  }
  return RgtaskxUB;
}
double compute_RgtaskLB(int i) {
  return max2double(compute_RgtaskxLB(i,'E'),compute_RgtaskxLB(i,'A'));
}
double compute_RgtaskUB(int i) {
  return max2double(compute_RgtaskxUB(i,'E'),compute_RgtaskxUB(i,'A'));
}

//   7. Functions to add lines to SMT file

void add_line_verbose100()        { fprintf(smt_instance_f,"(set-option :verbose 100)\n"); }
void add_line_maximize(char* obj) { fprintf(smt_instance_f,"(maximize %s)\n",obj); }
void add_line_checksat()          { fprintf(smt_instance_f,"(check-sat)\n"); }
void add_line_getobjectives()     { fprintf(smt_instance_f,"(get-objectives)\n"); }
void add_line_getmodel()          { fprintf(smt_instance_f,"(get-model)\n"); }
void createcommandstring_for_running_smt_solver() { fprintf(smt_cmd_f,"z3 %s > %s",smt_instance_fn,smt_result_fn); }

//   8. Functions for getting variable names

// vn means variable name.
void get_vnint1temp_for_rbfs_for_eq8(int i,int q,int j,char y,char x,int b,char* s) { sprintf(s,"int1temp_for_rbfs_for_eq8_%d_%d_%dUgU%cU%cU%d",i,q,j,y,x,b); }
void get_vnint2temp_for_rbfs_for_eq8(int i,int q,int j,char y,char x,int b,char* s) { sprintf(s,"int2temp_for_rbfs_for_eq8_%d_%d_%dUgU%cU%cU%d",i,q,j,y,x,b); }
void get_vnint1_for_rbfs_for_eq8(int i,int q,int j,char y,char x,int b,char* s) { sprintf(s,"int1_for_rbfs_for_eq8_%d_%d_%dUgU%cU%cU%d",i,q,j,y,x,b); }
void get_vnint2_for_rbfs_for_eq8(int i,int q,int j,char y,char x,int b,char* s) { sprintf(s,"int2_for_rbfs_for_eq8_%d_%d_%dUgU%cU%cU%d",i,q,j,y,x,b); }
void get_vnrbf_for_eq8(int i,int q,int j,char y,char x,int b,char* s) { sprintf(s,"rbf_for_eq8_%d_%d_%dUgU%cU%cU%d",i,q,j,y,x,b); }
void get_vnrbf_for_eq8phaseindicator(int i,int q,int j,char x,char* s) { sprintf(s,"rbf_for_eq8phaseindicator_%d_%d_%dUgU%c",i,q,j,x); }
void get_vnrbf_for_eq8phasevalue(int i,int q,int j,char x,char* s) { sprintf(s,"rbf_for_eq8phasevalue_%d_%d_%dUgU%c",i,q,j,x); }
void get_vnw(int i,int q,char x,char* s) { sprintf(s,"w_%d_%dUgU%c",i,q,x); }

//   9. Functions for adding variables

void add_variables_E() {
  int j;
  for (j=1;j<=ntasks;j++) {
    fprintf(smt_instance_f,"(declare-const E_%d Real)\n",j);
  }
}
void add_variables_integers_for_rbfs_for_eq8(int i,int q,char x) {
  int j; int iterator; char y; int b;
  char vn1temp[200]; char vn2temp[200];
  char vn1[200]; char vn2[200];
  for (j=1;j<=ntasks;j++) {
    for (iterator=0;iterator<=1;iterator++) {
      y = get_y_from_integer(iterator);
      for (b=0;b<=1;b++) {
        get_vnint1temp_for_rbfs_for_eq8(i,q,j,y,x,b,vn1temp);
        fprintf(smt_instance_f,"(declare-const %s Int)\n",vn1temp);
        get_vnint2temp_for_rbfs_for_eq8(i,q,j,y,x,b,vn2temp);
        fprintf(smt_instance_f,"(declare-const %s Int)\n",vn2temp);
        get_vnint1_for_rbfs_for_eq8(i,q,j,y,x,b,vn1);
        fprintf(smt_instance_f,"(declare-const %s Int)\n",vn1);
        get_vnint2_for_rbfs_for_eq8(i,q,j,y,x,b,vn2);
        fprintf(smt_instance_f,"(declare-const %s Int)\n",vn2);
      }
    }
  }
}
void add_variables_rbfs_for_eq8(int i,int q,char x) {
  int j; int iterator; char y; int b;
  char vn[200];
  for (j=1;j<=ntasks;j++) {
    for (iterator=0;iterator<=1;iterator++) {
      y = get_y_from_integer(iterator);
      for (b=0;b<=1;b++) {
        get_vnrbf_for_eq8(i,q,j,y,x,b,vn);
        fprintf(smt_instance_f,"(declare-const %s Real)\n",vn);
      }
    }
  }
}
void add_variable_rbf_for_eq8phaseindicator_and_value(int i,int q,char x) {
  int j;
  char vn[200];
  for (j=1;j<=ntasks;j++) {
    get_vnrbf_for_eq8phaseindicator(i,q,j,x,vn);
    fprintf(smt_instance_f,"(declare-const %s Bool)\n",vn);
    get_vnrbf_for_eq8phasevalue(i,q,j,x,vn);
    fprintf(smt_instance_f,"(declare-const %s Real)\n",vn);
  }
}
void add_variable_w(int i,int q,char x) {
  char vn[200];
  get_vnw(i,q,x,vn);
  fprintf(smt_instance_f,"(declare-const %s Real)\n",vn);
}

//  10. Functions for adding constraints

void add_constraints_Einrange() {
  int j;
  for (j=1;j<=ntasks;j++) {
    fprintf(smt_instance_f,"(assert ( <= %lf E_%d ))\n",LBE[j],j);    
    fprintf(smt_instance_f,"(assert ( <= E_%d %lf ))\n",j,UBE[j]);
  }
}
void add_constraint_ceiling_of_division(char* numeratorstr,char* denominatorstr,char* resultstr) {
  fprintf(smt_instance_f,"(assert ( < ( * ( - %s 1 ) %s ) %s ))\n",resultstr,denominatorstr, numeratorstr);
  fprintf(smt_instance_f,"(assert ( >= ( * %s %s ) %s ))\n",resultstr,denominatorstr,numeratorstr);
}
void add_constraint_non_negative_based_on_other(char* s_dest,char* s_source) {
  fprintf(smt_instance_f,"(assert ( => ( <  %s 0) ( = %s 0 ) ))\n",s_source,s_dest);
  fprintf(smt_instance_f,"(assert ( => ( >= %s 0) ( = %s %s) ))\n",s_source,s_dest,s_source);
}
void add_constraints_integers_for_rbfs_for_eq8(int i,int q,char x) {
  int j; int iterator; char y; int b;
  char vnw[200];
  char vn1temp[200]; char vn2temp[200];
  char vn1[200]; char vn2[200];
  char wstring[200];
  char integerstring1temp[200];
  char integerstring2temp[200];
  char integerstring1[200]; char numeratorintegerstring1[200]; char denominatorintegerstring1[200];
  char integerstring2[200]; char numeratorintegerstring2[200]; char denominatorintegerstring2[200];
  get_vnw(i,q,x,vnw);
  sprintf(wstring, "%s ",vnw);
  for (j=1;j<=ntasks;j++) {
    for (iterator=0;iterator<=1;iterator++) {
      y = get_y_from_integer(iterator);
      for (b=0;b<=1;b++) {
        get_vnint1temp_for_rbfs_for_eq8(i,q,j,y,x,b,vn1temp);
        get_vnint2temp_for_rbfs_for_eq8(i,q,j,y,x,b,vn2temp);
        sprintf(integerstring1temp,"%s ",vn1temp);
        sprintf(integerstring2temp,"%s ",vn2temp);
        get_vnint1_for_rbfs_for_eq8(i,q,j,y,x,b,vn1);
        get_vnint2_for_rbfs_for_eq8(i,q,j,y,x,b,vn2);
        sprintf(integerstring1,"%s ",vn1);
        sprintf(integerstring2,"%s ",vn2);
        if (y=='E') {
          sprintf(numeratorintegerstring1,   "( - %s ( - %lf E_%d ) ) ",wstring,tasks[j].T,j);
          sprintf(denominatorintegerstring1, "%lf ",tasks[j].T);  
          add_constraint_ceiling_of_division(numeratorintegerstring1,denominatorintegerstring1,integerstring1temp);
          add_constraint_non_negative_based_on_other(integerstring1,integerstring1temp);
          sprintf(numeratorintegerstring2,   "%s ",wstring);
          sprintf(denominatorintegerstring2, "%lf ",tasks[j].T);  
          add_constraint_ceiling_of_division(numeratorintegerstring2,denominatorintegerstring2,integerstring2);
        } else if (y=='A') {
          sprintf(numeratorintegerstring1,   "%s ",wstring);
          sprintf(denominatorintegerstring1, "%lf ",tasks[j].T);  
          add_constraint_ceiling_of_division(numeratorintegerstring1,denominatorintegerstring1,integerstring1);
          sprintf(numeratorintegerstring2,   "( - %s E_%d ) ",wstring,j);
          sprintf(denominatorintegerstring2, "%lf ",tasks[j].T);  
          add_constraint_ceiling_of_division(numeratorintegerstring2,denominatorintegerstring2,integerstring2temp);
          add_constraint_non_negative_based_on_other(integerstring2,integerstring2temp);
        } else {
            printf("Error in add_constraints_integers_for_rbfs_for_eq8\n"); exit(-1);
        }
      }
    }
  }
}
void add_constraints_rbfs_for_eq8(int i,int q,char x) {
  int j; int iterator; char y; int b;
  char vn[200];
  char vn1[200]; char vn2[200];
  char rbfstring[200];  
  char integerstring1[200];
  char integerstring2[200];
  for (j=1;j<=ntasks;j++) {
    for (iterator=0;iterator<=1;iterator++) {
      y = get_y_from_integer(iterator);
      for (b=0;b<=1;b++) {
        get_vnrbf_for_eq8(i,q,j,y,x,b,vn);
        sprintf(rbfstring, "%s ",vn);
        get_vnint1_for_rbfs_for_eq8(i,q,j,y,x,b,vn1);
        get_vnint2_for_rbfs_for_eq8(i,q,j,y,x,b,vn2);
        sprintf(integerstring1,"%s ",vn1);
        sprintf(integerstring2,"%s ",vn2);
        if        (b==1) {
          fprintf(smt_instance_f,"(assert ( = %s ( + ( * %s %lf ) ( * %s %lf ) ) ))\n",rbfstring,integerstring1,tasks[j].C,integerstring2,tasks[j].kC);
        } else if (b==0) {
          fprintf(smt_instance_f,"(assert ( = %s                  ( * %s %lf )   ))\n",rbfstring,             
                                        integerstring2,tasks[j].kC);       
        } else {
          printf("Error in add_constraints_rbfs_for_eq8\n"); exit(-1);        
        }
      }
    }
  }
}

void add_constraints_rbf_for_eq8phaseindicator_and_value(int i,int q,char x) {
  int j;
  char vnindicator[200];
  char vnvalue[200];
  char vnrbf1[200];
  char vnrbf2[200];
  for (j=1;j<=ntasks;j++) {
    get_vnrbf_for_eq8phaseindicator(i,q,j,x,vnindicator);
    get_vnrbf_for_eq8phasevalue(i,q,j,x,vnvalue);
    get_vnrbf_for_eq8(i,q,j,'E',x,1,vnrbf1);
    get_vnrbf_for_eq8(i,q,j,'A',x,1,vnrbf2);
    fprintf(smt_instance_f,"(assert ( => ( not %s ) ( >= %s %s) ))\n",vnindicator,vnrbf1,vnrbf2);
    fprintf(smt_instance_f,"(assert ( =>       %s   ( >= %s %s) ))\n",vnindicator,vnrbf2,vnrbf1);
    fprintf(smt_instance_f,"(assert ( => ( not %s ) ( =  %s %s) ))\n",vnindicator,vnvalue,vnrbf1);
    fprintf(smt_instance_f,"(assert ( =>       %s   ( =  %s %s) ))\n",vnindicator,vnvalue,vnrbf2);
  }
}
void add_constraints_eq8(int i,int q,char x) {
  int j; char y; int b;
  char vnw[200];
  char vnrbf[200];
  char wstring[200];
  char vnvalue[200];
  if (ntasks<=0) { printf("Error in add_constraints_eq8\n"); exit(-1); }
  get_vnw(i,q,x,vnw);
  sprintf(wstring, "%s ",vnw);
  fprintf(smt_instance_f,"(assert ( = %s ( + ",wstring);
  for (j=1;j<=ntasks;j++) {        
    if (tasks[j].prio<tasks[i].prio) {
      y = 'E'; b = 0;
      get_vnrbf_for_eq8(i,q,j,y,x,b,vnrbf);
      fprintf(smt_instance_f,"%s ",vnrbf);   
    }
  }
  if        (x=='E') {
    fprintf(smt_instance_f,"%lf ",q*tasks[i].C+q*tasks[i].kC);
  } else if (x=='A') {
    fprintf(smt_instance_f,"%lf ",q*tasks[i].C+(q-1)*tasks[i].kC);     
  }
  for (j=1;j<=ntasks;j++) {        
    if (tasks[j].prio>tasks[i].prio) {
      get_vnrbf_for_eq8phasevalue(i,q,j,x,vnvalue);
      fprintf(smt_instance_f, "%s ",vnvalue);
    }
  }
  fprintf(smt_instance_f, "0 ) )) \n"); // We want to make sure that the summation over terms has at least
                                                       // two terms. if ntasks=1 and j=1, then we would end up with just one
                                                       // term. But since we have 0 here, we end up with two terms for this
                                                       // case as well.
}
void add_constraints_line26(int i,int q,char x) {
  char vn[200];
  char wstring[200];
  get_vnw(i,q,x,vn);
  sprintf(wstring, "%s ",vn);
  fprintf(smt_instance_f,"(assert ( <= ( - %s ",wstring);
  if        (x=='E') {
    fprintf(smt_instance_f, "( + %lf (- %lf E_%d) ) ", (q-1)*tasks[i].T,tasks[i].T,i);
  } else if (x=='A'){
    fprintf(smt_instance_f, "%lf ",(q-1)*tasks[i].T);
  } else {
    printf("Error in add_constraints_line26\n"); exit(-1);
  }
  fprintf(smt_instance_f, ") E_%d ))\n",i);
}

//  11. Functions for retrieving information from result file from SMT solving

char* get_first_non_blank_char(char* str) {
  char* p;
  p = str;
  while ((*p)!=0) {
    if (char_is_blank(*p)) {
      p++;
    } else {
      return p;    
    }
  }
  return NULL;    
}

// It is assumed that for the string str, all leading and trailing blanks have been deleted
// and also that all outer parentheses have been deleted.
double retrieve_real_number_from_string(char* str) {
  int n_elements_read; double v; double num; double denom; char* p;
  n_elements_read = sscanf(str,"%lf",&v);
  if (n_elements_read==1) {
    return v;  
  } else {
    p = get_first_non_blank_char(str);
    if (p!=NULL) {
      if ((*p)=='-')  {
        return (-1.0)*retrieve_real_number_from_string(p+1);   
      } else {
        if ((*p)=='/')  {
          n_elements_read = sscanf(p+1,"%lf %lf",&num,&denom);        
          if (n_elements_read==2) {
            return num/denom;
          } else {
            printf("Error in retrieve_real_number_from_string. Expected to get two reals. str=%s p=%s\n",str,p); exit(-1);
          }
        } else {
          printf("Error in retrieve_real_number_from_string. Expected to get division. str=%s p=%s\n",str,p); exit(-1);
        }
      }
    } else {
      printf("Error in retrieve_real_number_from_string. Expected to get a non-blank character. str=%s\n",str); exit(-1);    
    }
  }
}

// This function assumes that smt_result_fn has data and is satisfiable.
double getvalueofrealvariable_from_results_file(char* varname) {
  char s[200]; char s2[200]; char* p; char* p2; double v;
  FILE* smt_result_f;
  smt_result_f = fopen(smt_result_fn,"r");
  if (smt_result_f==NULL) { printf("In getvalueofrealvariable_from_results_file. Error when opening file. fn=%s\n",smt_instance_fn); exit(0); }
  p  = fgets(s,200,smt_result_f);
  if (p==NULL) { printf("In getvalueofrealvariable_from_results_file. Error when reading 1st line. fn=%s\n",smt_instance_fn); exit(0); }
  p2 = fgets(s2,200,smt_result_f);
  if (p2==NULL) { printf("In getvalueofrealvariable_from_results_file. Error when reading 2nd line. fn=%s\n",smt_instance_fn); exit(0); }
  while (p2!=NULL) {
    if (strstr(s,varname)!=NULL) {
      trim_spaces_and_parentheses(s2);
      v = retrieve_real_number_from_string(s2);
      return v;
    }  
    strcpy(s,s2);
    p2 = fgets(s2,200,smt_result_f);
  }
  fclose(smt_result_f);
}

//  12. Functions for bounds on E

void print_LBE_and_UBE() {
  int j;
  for (j=1;j<=ntasks;j++) {
    printf("LBE[%d] = %lf UBE[%d] = %lf\n",j,LBE[j],j,UBE[j]);
  }
}
void initialize_all_LBE_and_UBE(double* RK) {
  int j;
  for (j=1;j<=ntasks;j++) {
    LBE[j] = tasks[j].C;
    UBE[j] = tasks[j].D-RK[j];
  }
}
int calc_potentially_feasible_from_LBE_and_UBE() {
  int j;
  for (j=1;j<=ntasks;j++) {
    if (UBE[j]<LBE[j]) {
      return 0;    
    }
  }
  return 1;    
}

int compute_all_LBRg_and_then_update_all_LBE() {
  int changed; int i; double LBRg_this_task;
  changed = 0;
  for (i=1;i<=ntasks;i++) {
    LBRg_this_task = compute_RgtaskLB(i);
    if (LBE[i]<LBRg_this_task) {
      LBE[i] = LBRg_this_task; changed = 1;
    }
  }
  return changed;
}

//  13. Functions for solving PR

void declare_all_variables_solve_PR_satisfaction(double* tUBUgUE,double* tUBUgUA) {
  int i; int q; int iterator; char x; double t;
  add_variables_E();
  for (i=1;i<=ntasks;i++) {
    for (iterator=0;iterator<=1;iterator++) {
      x = get_x_from_integer(iterator);
      if      (x=='E') { t = tUBUgUE[i]; }
      else if (x=='A') { t = tUBUgUA[i]; }
      else             { printf("Error in declare_all_variables_solve_PR_satisfaction.\n"); exit(-1); }
      for (q=1;q<=ceil(t/tasks[i].T);q++) {
        add_variables_integers_for_rbfs_for_eq8(i,q,x);
        add_variables_rbfs_for_eq8(i,q,x);
        add_variable_rbf_for_eq8phaseindicator_and_value(i,q,x);
        add_variable_w(i,q,x);
      }
    }
  }
}

void define_all_constraints_solve_PR_satisfaction(double* tUBUgUE,double* tUBUgUA) {
  int i; int q; int iterator; char x; double t;
  add_constraints_Einrange();
  for (i=1;i<=ntasks;i++) {
    for (iterator=0;iterator<=1;iterator++) {
      x = get_x_from_integer(iterator);
      if      (x=='E') { t = tUBUgUE[i]; }
      else if (x=='A') { t = tUBUgUA[i]; }
      else             { printf("Error in define_all_constraints_solve_PR_satisfaction. First loop.\n"); exit(-1); }
      for (q=1;q<=ceil(t/tasks[i].T);q++) {
        add_constraints_integers_for_rbfs_for_eq8(i,q,x);
        add_constraints_rbfs_for_eq8(i,q,x);
        add_constraints_rbf_for_eq8phaseindicator_and_value(i,q,x);
        add_constraints_eq8(i,q,x);
      }
    }
  }
  for (i=1;i<=ntasks;i++) {
    for (iterator=0;iterator<=1;iterator++) {
      x = get_x_from_integer(iterator);
      if      (x=='E') { t = tUBUgUE[i]; }
      else if (x=='A') { t = tUBUgUA[i]; }
      else             { printf("Error in define_all_constraints_solve_PR_satisfaction. Second loop.\n"); exit(-1); }
      for (q=1;q<=ceil(t/tasks[i].T);q++) {
        add_constraints_line26(i,q,x);
      }
    }
  }
}

int solve_PR_satisfaction(double* tUBUgUE,double* tUBUgUA,double* RK,double* E_values) {
  char temps[200]; int j;
  counterforSMTsolver = counterforSMTsolver + 1;
  sprintf(smt_instance_fn,"smt_instance_fn_%d.smt2",counterforSMTsolver);
  sprintf(smt_cmd_fn,"./smt_cmd_%d.cmd",counterforSMTsolver);
  sprintf(smt_result_fn,"smt_result_%d.txt",counterforSMTsolver);
  smt_instance_f = fopen(smt_instance_fn,"w");
  if (smt_instance_f==NULL) { printf("In solve_PR_satisfaction. Error when opening file. fn=%s\n",smt_instance_fn); exit(0); }
  smt_cmd_f = fopen(smt_cmd_fn,"w");
  if (smt_cmd_f==NULL) { printf("In solve_PR_satisfaction. Error when opening file. fn=%s\n",smt_cmd_fn); exit(0); }
  // if (VERBOSE) { add_line_verbose100(); }
  declare_all_variables_solve_PR_satisfaction(tUBUgUE,tUBUgUA);
  define_all_constraints_solve_PR_satisfaction(tUBUgUE,tUBUgUA);
  add_line_checksat();
  add_line_getmodel();
  createcommandstring_for_running_smt_solver();
  fclose(smt_instance_f);
  fclose(smt_cmd_f);
  setuserexecutablemode(smt_cmd_fn);
  system(smt_cmd_fn);  
  if (is_the_smt_result_satisfiable()) {
    for (j=1;j<=ntasks;j++) {
      sprintf(temps,"E_%d",j); 
      E_values[j] = getvalueofrealvariable_from_results_file(temps);
    }
    return 1;
  } else {
    return 0;
  }
}

//  14. Function for optimal configuration for assigning only E

void compute_tUBUgUE_and_tUBUgUA(double* tUBUgUE, double* tUBUgUA) {
  int i; int iterator; char x;
  for (i=1;i<=ntasks;i++) {
    for (iterator=0;iterator<=1;iterator++) {
      x = get_x_from_integer(iterator);
      if        (x=='E') {
        tUBUgUE[i] = compute_tUB(i,x);          
        if (VERBOSE) { printf("tUBUgUE[%d] = %lf\n",i,tUBUgUE[i]); }
      } else if (x=='A') {
        tUBUgUA[i] = compute_tUB(i,x);
        if (VERBOSE) { printf("tUBUgUA[%d] = %lf\n",i,tUBUgUA[i]); }
      } else {
        printf("Error in perform_optimal_configuration_assign_only_E.\n"); exit(-1);        
      }
    }
  }
}
int perform_optimal_configuration_assign_only_E() { 
  int ret_val; int flag; int changed; int i;
  double* RK; double* tUBUgUE; double* tUBUgUA; double* computed_E_values;
  if (!priorities_are_unique()) { printf("Failure in perform_optimal_configuration_assign_only_E. We assume that priorities are unique but this is not the case."); exit(-1); }  
  RK = malloc(sizeof(double)*(ntasks+1)); if (RK==NULL) { printf("Memory allocation failed in perform_optimal_configuration_assign_only_E."); exit(-1); }
  tUBUgUE = malloc(sizeof(double)*(ntasks+1)); if (tUBUgUE==NULL) { printf("Memory allocation failed in perform_optimal_configuration_assign_only_E."); exit(-1); }
  tUBUgUA = malloc(sizeof(double)*(ntasks+1)); if (tUBUgUA==NULL) { printf("Memory allocation failed in perform_optimal_configuration_assign_only_E."); exit(-1); }
  computed_E_values = malloc(sizeof(double)*(ntasks+1)); if (computed_E_values==NULL) { printf("Memory allocation failed in perform_optimal_configuration_assign_only_E."); exit(-1); }
  if (VERBOSE) { printf("getutil() = %lf\n",getutil()); fflush(stdout); }
  if (getutil()<=1.0-UTIL_EPSILON) {
    if (VERBOSE) { printf("Starting analyzing hypertasks\n"); fflush(stdout); }
    for (i=1;i<=ntasks;i++) {
      RK[i] = compute_hypertask_responsetime(i);
    }
    if (VERBOSE) { printf("Done analyzing hypertasks\n"); fflush(stdout); }    
    if (hypertask_response_times_respect_deadlines(RK)) {
      if (VERBOSE) { printf("Initializing bounds on E values\n"); fflush(stdout); }
      initialize_all_LBE_and_UBE(RK);
      if (VERBOSE) { printf("Updating bounds on E values\n"); fflush(stdout); }
      changed = compute_all_LBRg_and_then_update_all_LBE();
      while ((changed) && (calc_potentially_feasible_from_LBE_and_UBE())) {
        changed = compute_all_LBRg_and_then_update_all_LBE();
      }    
      if (VERBOSE) { printf("Done computing bounds on E values\n"); fflush(stdout); }
      if (VERBOSE) { print_LBE_and_UBE(); }
      if (calc_potentially_feasible_from_LBE_and_UBE()) {
        if (VERBOSE) { printf("Starting computing upper bound on busy period\n"); fflush(stdout); }
        compute_tUBUgUE_and_tUBUgUA(tUBUgUE,tUBUgUA);
        if (VERBOSE) { printf("Done computing upper bound on busy period\n"); fflush(stdout); }
        if (VERBOSE) { printf("Starting computing E values\n"); fflush(stdout); }
        flag = solve_PR_satisfaction(tUBUgUE,tUBUgUA,RK,computed_E_values);
        if (VERBOSE) { printf("Done computing E values\n"); fflush(stdout); }
        if (flag) {
          for (i=1;i<=ntasks;i++) {
            tasks[i].E = computed_E_values[i]; 
          }
          printE_to_fn(outputfilename);
          ret_val = 1;
        } else {
          ret_val = 0;
        }
      } else {
        ret_val = 0;
      }
    } else {
      ret_val = 0;   
    }
  } else {
    ret_val = 0;
  }
  free(computed_E_values);
  free(tUBUgUE);
  free(tUBUgUA);
  free(RK);
  return ret_val;
}

//  15. Functions for generating random taskset

void adjustCkC_to_each_targetutil(double targetutil) {
  double scalefactor; int j;
  scalefactor = targetutil/getutil();
  for (j=1;j<=ntasks;j++) {
    tasks[j].C  = tasks[j].C  * scalefactor;
    tasks[j].kC = tasks[j].kC * scalefactor;
  }
}
void sort_tasks_so_that_smallest_deadline_has_lowest_index() {
  int j; int jprime;
  struct task temp;
  for (j=1;j<=ntasks-1;j++) {
    for (jprime=j+1;jprime<=ntasks;jprime++) {
      if (tasks[j].D>tasks[jprime].D) {
        temp = tasks[j];
        tasks[j] = tasks[jprime];
        tasks[jprime] = temp;
      }
    }
  }
}
void generate_taskset_randomly(int nt,int TUTILINDEX,int TMAXEXP) {
  int j; double CTOT;
  usemethod = 1;
  ntasks = nt;
  for (j=1;j<=ntasks;j++) {
    tasks[j].T  = mypow(2.0,my_double_random(0.0,log2(TMAXEXP)));  
    tasks[j].D  = 0.8*tasks[j].T;
    tasks[j].E  = 0.8*tasks[j].D;
    tasks[j].C  = 0.9*tasks[j].D;
    tasks[j].kC = 0.1*tasks[j].D;
  }
  adjustCkC_to_each_targetutil(TUTILINDEX/10.0); // we permit generating tasksets with utilization that is exactly 1;
                                                 // in schedulability analysis of such tasksets, we will declare
                                                 // 'can't guarantee' early.
  sort_tasks_so_that_smallest_deadline_has_lowest_index();
  for (j=1;j<=ntasks;j++) {
    tasks[j].id = j;
    tasks[j].prio = ntasks+1-j;
  }
}
void get_file_name_for_taskset(int nt,int TUTILINDEX,int TMAXEXP,char* fn_for_taskset) {
  sprintf(fn_for_taskset,"taskset_%d_%d_%d.txt",nt,TUTILINDEX,TMAXEXP);
}
void get_file_name_for_taskset_outputfilename(int nt,int TUTILINDEX,int TMAXEXP,char* outp_fn_for_taskset) {
  sprintf(outp_fn_for_taskset,"output_taskset_%d_%d_%d.txt",nt,TUTILINDEX,TMAXEXP);
}
void do_metathing() {
  struct timespec mybegin;
  struct timespec myend;
  int status_clock_gettime_begin;
  int status_clock_gettime_end;
  int nt; int TUTILINDEX; int TMAXEXP; int DTRATIOINDEX;
  int flag; double timerequired;
  char fn_for_taskset[2000];
  char outp_fn_for_taskset[2000];
  FILE* myfile;
  myfile = fopen("overallstatistics.txt","w");
  for (nt=2;nt<=5;nt++) {
    for (TUTILINDEX=1;TUTILINDEX<=11;TUTILINDEX++) {
      for (TMAXEXP=1;TMAXEXP<=128;TMAXEXP*=2) {
        generate_taskset_randomly(nt,TUTILINDEX,TMAXEXP);
        get_file_name_for_taskset(nt,TUTILINDEX,TMAXEXP,fn_for_taskset);
        writetasksettofile(fn_for_taskset);
        status_clock_gettime_begin = clock_gettime(CLOCK_MONOTONIC_RAW, &mybegin);
        flag = perform_optimal_configuration_assign_only_E();
        status_clock_gettime_end = clock_gettime(CLOCK_MONOTONIC_RAW, &myend);
        timerequired = diff_s_timespec(&myend,&mybegin);
        fprintf(myfile,"%d %d %d  %d %lf\n",nt,TUTILINDEX,TMAXEXP,flag,timerequired );fflush(myfile);
        get_file_name_for_taskset_outputfilename(nt,TUTILINDEX,TMAXEXP,outp_fn_for_taskset);
        printE_to_fn(outp_fn_for_taskset);
      }      
    }    
  }
  fclose(myfile);
}

//  16. Main

int main( int argc, char *argv[], char *envp[]) {
  int flag;
  initialize_random_generator();
  parsearguments(argc, argv, envp);
  if (metathing) {
    do_metathing();
    exit(0);
  }
  readtasksetfromfile( inputfilename );
  flag = perform_optimal_configuration_assign_only_E();
  if (flag) {
    printf("Taskset is schedulable.\n");
    return 1;
  } else {
    printf("Taskset is unschedulable.\n");
    return 0;
  }
}

//  17. Test results
// 
// ba@ba-XPS-8900:~/mtoe$ cat system.txt system2.txt system3.txt 
// 1
// 2
// 1 2
//   20 20 10 8 3
// 2 1
//   25 25 12 7 2
// 
// 1
// 2
// 1 2
//   20 20 10 9 3
// 2 1
//   25 25 12 7 2
// 
// 1
// 2
// 1 2
//   20 20 10 8 3
// 2 1
//   25 25 12 8 2
// 
// ba@ba-XPS-8900:~/mtoe/fixing$ ./mtoe -i system.txt -o system_results_from_analysis.txt
// Taskset is schedulable.
// ba@ba-XPS-8900:~/mtoe/fixing$ ./mtoe -i system2.txt -o system2_results_from_analysis.txt
// Taskset is unschedulable.
// ba@ba-XPS-8900:~/mtoe/fixing$ ./mtoe -i system3.txt -o system3_results_from_analysis.txt
// Taskset is unschedulable.
// ba@ba-XPS-8900:~/mtoe/fixing$ cat system_results_from_analysis.txt
// E_1 = 10.000000
// E_2 = 20.000000