24.09., 9:00 - 11:00: Due to updates GitLab will be unavailable for some minutes between 09:00 and 11:00.

Commit 55128663 authored by Pushpak Jagtap's avatar Pushpak Jagtap Committed by GitHub

Add files via upload

parent 1f981577
#
# compiler
#
#CC = g++
CC = clang++
CXXFLAGS = -Wall -Wextra -std=c++11 -O3 -g -O0
#
# scots
#
SCOTSROOT = ../../..
SCOTSINC = -I$(SCOTSROOT)/bdd -I$(SCOTSROOT)/utils
#
# cudd
#
CUDDPATH = /opt/local
CUDDINC = -I$(CUDDPATH)/include
CUDDLIBS = -lcudd
CUDDLPATH = -L$(CUDDPATH)/lib
TARGET = thermalmodel10R
all: $(TARGET)
%.o:%.cc
$(CC) -c $(CXXFLAGS) $(CUDDINC) $(SCOTSINC) $< -o $@
$(TARGET): $(TARGET).o
$(CC) $(CXXFLAGS) -o $(TARGET) $(TARGET).o $(CUDDLPATH) $(CUDDLIBS)
clean:
rm ./$(TARGET) ./$(TARGET).o
#include <array>
#include <iostream>
#include<cmath>
#include "cuddObj.hh"
#include "../../src/TicToc.hh"
#include "../../src/fixedPointMode.hh"
#include "../../src/getAbstraction.hh"
#include "../../src/SymbolicSetSpace.hh"
#include "../../src/abstractionMode.hh"
using namespace std;
const int sDIM = 10; /* System dimension */
const int iDIM = 2; /* Input dimension */
const double T = 25; /* Sampling time */
size_t P = 3; /* Cardinality of (quantized) input set */
size_t N = 12; /* Temporal Horizon */
typedef std::array<double,sDIM> state_type;
/* forward declaration of the ode solver */
template<class F>
void ode_solver(F rhs, state_type &x, size_t nint, double h);
/* ODE */
auto system_post = [](state_type &x, int u) -> void {
/* Mapping abstract input to system inputs */
double us[iDIM]={0};
if(u==0)
{us[0]=0;us[1]=0;}
if(u==1)
{us[0]=0;us[1]=1;}
if(u==2)
{us[0]=1;us[1]=0;}
/* ode describing thermal model*/
auto rhs=[us](state_type &xx, const state_type &x) -> void {
const double a=0.05;
const double ae2=0.005;
const double ae5=0.005;
const double ae=0.0033;
const double ah=0.0036;
const double te=12;
const double th=100;
xx[0]=(-a-ae)*x[0]+a*x[1]+ae*te;
xx[1]=(-4*a-ae2-ah*us[0])*x[1]+a*x[0]+a*x[6]+a*x[8]+a*x[2]+ae2*te+ah*th*us[0];
xx[2]=(-2*a-ae)*x[2]+a*x[1]+a*x[3]+ae*te;
xx[3]=(-2*a-ae)*x[3]+a*x[2]+a*x[4]+ae*te;
xx[4]=(-4*a-ae5-ah*us[1])*x[4]+a*x[3]+a*x[7]+a*x[5]+a*x[9]+ae5*te+ah*th*us[1];
xx[5]=(-a-ae)*x[5]+a*x[4]+ae*te;
xx[6]=(-a-ae)*x[6]+a*x[1]+ae*te;
xx[7]=(-a-ae)*x[7]+a*x[4]+ae*te;
xx[8]=(-a-ae)*x[8]+a*x[1]+ae*te;
xx[9]=(-a-ae)*x[9]+a*x[4]+ae*te;
};
size_t nint = 5; /* no. of time step for ode solving */
double h=T/nint; /* time step for ode solving (T is an sampling time) */
ode_solver(rhs,x,nint,h); /* Runga Kutte solver */
};
/* defining constrains for the controller
ul[i] : upper limit for the temperature in ith dimension
ll[i] : lower limit for the temperature in ith dimension */
auto constrain = [](state_type y) -> bool {
double ul=21.8;
double ll=18;
bool s = true;
for(int j = 0; j < sDIM; j++)
{
if( y[j] >= ul || y[j] <= ll )
{
s = false;
break;
}
}
return s;
};
int main(){
/*to measure time */
TicToc tt;
/* CUDD Manager */
Cudd ddmgr;
/* source state of the system */
state_type xs;
xs[0]=17;xs[1]=17;xs[2]=17;xs[3]=17;xs[4]=17;xs[5]=17;xs[6]=17;xs[7]=17;xs[8]=17;xs[9]=17;
/* initial state of the system */
state_type x0;
x0[0]=19;x0[1]=19;x0[2]=19;x0[3]=19.8;x0[4]=20.8;x0[5]=19.8;x0[6]=19;x0[7]=19.8;x0[8]=19;x0[9]=19.8;
/* approximation on abstraction */
double epsilon = 0.1;
/* number of time samples the controller should run */
size_t t = 20;
/******************************************************
Symbolic model construction
******************************************************/
/* defining Symbolic Set for the transition relation */
SymbolicSetSpace ss(ddmgr,P,N);
tt.tic();
ss.addAbsStates();
/* defining abstraction class */
getAbstraction<state_type> ab(&ss);
tt.toc();
std::cout<<"Number of Elements in the Transition Relation : "<<ab.getSize()<<std::endl;
/* Computing the constrain set */
std::cout<<"Computing the constrain set ... ";
tt.tic();
BDD set = ab.getConstrainSet(system_post,constrain,xs);
tt.toc();
if(set == ddmgr.bddZero()){
std::cout << "Set is empty !!";
return 0;
}
/******************************************************
synthesizing the controller
******************************************************/
std::cout<<"Computing the controller ... ";
fixedPointMode fp(&ab);
BDD C;
tt.tic();
/* controller for safety specification */
C = fp.safe(set);
tt.toc();
if(C == ddmgr.bddZero()){
std::cout << "No controller is found !!";
return 0;
}
/* find out mode sequence which is very close to the output */
std::cout<<"Computed initial state is...\n";
/* finding inital stae in abstraction for safety specification */
BDD w0 = ab.getAbsState(system_post,set,x0,xs,epsilon,sDIM);
/* open-loop simulation */
ab.openLoopSim(C,w0,system_post,x0,sDIM,t);
return 0;
}
template<class F>
void ode_solver(F rhs, state_type &x, size_t nint, double h) {
/* runge kutte order 4 */
state_type k[4];
state_type tmp;
for(size_t t=0; t<nint; t++) {
rhs(k[0],x);
for(size_t i=0;i<sDIM;i++)
tmp[i]=x[i]+h/2*k[0][i];
rhs(k[1],tmp);
for(size_t i=0;i<sDIM;i++)
tmp[i]=x[i]+h/2*k[1][i];
rhs(k[2],tmp);
for(size_t i=0;i<sDIM;i++)
tmp[i]=x[i]+h*k[2][i];
rhs(k[3],tmp);
for(size_t i=0; i<sDIM; i++)
x[i] = x[i] + (h/6)*(k[0][i] + 2*k[1][i] + 2*k[2][i] + k[3][i]);
}
}
clc;
clear all;
sDIM=10;
a=5e-2;
ae2=5e-3;ae5=5e-3;
ae1=3.3e-3;ae3=3.3e-3;ae4=3.3e-3;ae6=3.3e-3;ae7=3.3e-3; ae8=3.3e-3; ae9=3.3e-3;ae10=3.3e-3;
ah=3.6e-3;
te=12;
th=100;
Dt=5; %time step for ode solving
X1=zeros(sDIM,1);
%control inputs
filename = 'controller.txt';
input=importdata(filename);
for i=1:length(input)
u(1,i)=bitand(floor(input(i)/2),1);
u(2,i)=bitand(input(i),1);
end
% initial conditions
X_init=[18.9 19 19.1 19.5 20.8 19.7 19 19.8 19 19.8];
%initialization
X=X_init;
time = 0;
k=0;
for j = 1:length(input)*5
k=ceil(j/5);
X1(1)=X(1)+Dt*((-a-ae1)*X(1)+a*X(2)+ae1*te);
X1(2)=X(2)+Dt*((-4*a-ae2-ah*u(1,k))*X(2)+a*X(1)+a*X(7)+a*X(9)+a*X(3)+ae2*te+ah*th*u(1,k));
X1(3)=X(3)+Dt*((-2*a-ae3)*X(3)+a*X(2)+a*X(4)+ae3*te);
X1(4)=X(4)+Dt*((-2*a-ae4)*X(4)+a*X(3)+a*X(5)+ae4*te);
X1(5)=X(5)+Dt*((-4*a-ae5-ah*u(2,k))*X(5)+a*X(4)+a*X(8)+a*X(6)+a*X(10)+ae5*te+ah*th*u(2,k));
X1(6)=X(6)+Dt*((-a-ae6)*X(6)+a*X(5)+ae6*te);
X1(7)=X(7)+Dt*((-a-ae7)*X(7)+a*X(2)+ae7*te);
X1(8)=X(8)+Dt*((-a-ae8)*X(8)+a*X(5)+ae8*te);
X1(9)=X(9)+Dt*((-a-ae9)*X(9)+a*X(2)+ae9*te);
X1(10)=X(10)+Dt*((-a-ae10)*X(10)+a*X(5)+ae10*te);
X=X1;
Y(1,j)=X(1);
Y(2,j)=X(2);
Y(3,j)=X(3);
Y(4,j)=X(4);
Y(5,j)=X(5);
Y(6,j)=X(6);
Y(7,j)=X(7);
Y(8,j)=X(8);
Y(9,j)=X(9);
Y(10,j)=X(10);
Time(1,j)=time;
time=time+Dt;
end
time1=0;
for l = 1:length(input)*500
m=ceil(l/500);
Uem(1,l)=u(1,m);
Uem(2,l)=u(2,m);
time1=time1+Dt/100;
Time1(1,l)=time1;
end
figure(1);
plot(Time,Y);
title('Output Response');
figure(2);
plot(Time1,Uem(1,:));
title('control input to heater 1');
figure(3);
plot(Time1,Uem(2,:));
title('control input to heater 2');
#
# compiler
#
#CC = g++
CC = clang++
CXXFLAGS = -Wall -Wextra -std=c++11 -O3 -g -O0
#
# scots
#
SCOTSROOT = ../../..
SCOTSINC = -I$(SCOTSROOT)/bdd -I$(SCOTSROOT)/utils
#
# cudd
#
CUDDPATH = /opt/local
CUDDINC = -I$(CUDDPATH)/include
CUDDLIBS = -lcudd
CUDDLPATH = -L$(CUDDPATH)/lib
TARGET = thermalmodel6R
all: $(TARGET)
%.o:%.cc
$(CC) -c $(CXXFLAGS) $(CUDDINC) $(SCOTSINC) $< -o $@
$(TARGET): $(TARGET).o
$(CC) $(CXXFLAGS) -o $(TARGET) $(TARGET).o $(CUDDLPATH) $(CUDDLIBS)
clean:
rm ./$(TARGET) ./$(TARGET).o
#include <array>
#include <iostream>
#include<cmath>
#include "cuddObj.hh"
#include "../../src/TicToc.hh"
#include "../../src/fixedPointMode.hh"
#include "../../src/getAbstraction.hh"
#include "../../src/SymbolicSetSpace.hh"
#include "../../src/abstractionMode.hh"
using namespace std;
const int sDIM = 6; /* System dimension */
const int iDIM = 2; /* Input dimension */
const double T = 25; /* Sampling time */
size_t P = 3; /* Cardinality of (quantized) input set */
size_t N = 12; /* Temporal Horizon */
typedef std::array<double,sDIM> state_type;
/* forward declaration of the ode solver */
template<class F>
void ode_solver(F rhs, state_type &x, size_t nint, double h);
/* ODE */
auto system_post = [](state_type &x, int u) -> void {
/* Defining system inputs to system*/
double ub[iDIM]={0};
if(u==0)
{ub[0]=0;ub[1]=0;}
if(u==1)
{ub[0]=0;ub[1]=1;}
if(u==2)
{ub[0]=1;ub[1]=0;}
/* ode describing thermal model*/
auto rhs=[ub](state_type &xx, const state_type &x) -> void {
const double a=0.05;
const double ae1=0.005;
const double ae4=0.005;
const double ae=0.0033;
const double ah=0.0036;
const double te=12;
const double th=100;
xx[0] = (-3*a-ae1-ah*ub[0])*x[0]+a*x[1]+a*x[2]+a*x[4]+ae1*te+ah*th*ub[0];
xx[1] = (-2*a-ae)*x[1]+a*x[0]+a*x[3]+ae*te;
xx[2] = (-2*a-ae)*x[2]+a*x[0]+a*x[3]+ae*te;
xx[3] = (-3*a-ae4-ah*ub[1])*x[3]+a*x[1]+a*x[2]+a*x[5]+ae4*te+ah*th*ub[1];
xx[4] = (-a-ae)*x[4]+a*x[0]+ae*te;
xx[5] = (-a-ae)*x[5]+a*x[3]+ae*te;
};
size_t nint = 5; /* no. of time step for ode_solving */
double h=T/nint;
ode_solver(rhs,x,nint,h);
};
/* defining constrains for the controller
ul[i] : upper limit of the ith dimension
ll[i] : lower limit of the ith dimension */
auto constrain = [](state_type y) -> bool {
double ul=21.0;
double ll=17.5;
bool s = true;
for(int j = 0; j < sDIM; j++)
{
if( y[j] >= ul || y[j] <= ll )
{
s = false;
break;
}
}
return s;
};
int main(){
/*to measure time */
TicToc tt;
/* CUDD Manager */
Cudd ddmgr;
/* source state of the system */
state_type xs;
xs[0]=17;xs[1]=17;xs[2]=17;xs[3]=17;xs[4]=17;xs[5]=17;
/* initial state of the system */
state_type x0;
x0[0]=16;x0[1]=16;x0[2]=16;x0[3]=16;x0[4]=16;x0[5]=16;
/* approximation on abstraction */
double epsilon = 0.2;
/* number of time samples the controller should run */
size_t t = 40;
/******************************************************
Symbolic model construction
******************************************************/
/* defining Symbolic Set for the transition relation */
SymbolicSetSpace ss(ddmgr,P,N);
tt.tic();
ss.addAbsStates();
/* defining abstraction class */
getAbstraction<state_type> ab(&ss);
tt.toc();
std::cout<<"Number of Elements in the Transition Relation : "<<ab.getSize()<<std::endl;
/* Computing the constrain set */
std::cout<<"Computing the constrain set ... ";
tt.tic();
BDD set = ab.getConstrainSet(system_post,constrain,xs);
tt.toc();
if(set == ddmgr.bddZero()){
std::cout << "Set is empty !!";
return 0;
}
/******************************************************
synthesizing the controller
******************************************************/
std::cout<<"Computing the controller ... ";
fixedPointMode fp(&ab);
BDD C;
tt.tic();
/* controller for reach and stay specification */
C = fp.reachStay(set);
tt.toc();
if(C == ddmgr.bddZero()){
std::cout << "No controller is found !!";
return 0;
}
/* find out mode sequence which is very close to the output */
std::cout<<"Computed the initial state is...\n";
/* finding inital stae in abstraction for safety specification */
BDD w0 = ab.getAbsState(system_post,x0,xs,epsilon,sDIM);
/* open-loop simulation */
ab.openLoopSim(C,w0,system_post,x0,sDIM,t);
return 0;
}
template<class F>
void ode_solver(F rhs, state_type &x, size_t nint, double h) {
/* runge kutte order 4 */
state_type k[4];
state_type tmp;
for(size_t t=0; t<nint; t++) {
rhs(k[0],x);
for(size_t i=0;i<sDIM;i++)
tmp[i]=x[i]+h/2*k[0][i];
rhs(k[1],tmp);
for(size_t i=0;i<sDIM;i++)
tmp[i]=x[i]+h/2*k[1][i];
rhs(k[2],tmp);
for(size_t i=0;i<sDIM;i++)
tmp[i]=x[i]+h*k[2][i];
rhs(k[3],tmp);
for(size_t i=0; i<sDIM; i++)
x[i] = x[i] + (h/6)*(k[0][i] + 2*k[1][i] + 2*k[2][i] + k[3][i]);
}
}
clc;
clear all;
sDIM=6;
a=5e-2;
ae1=5e-3;ae4=5e-3;
ae2=3.3e-3;ae3=3.3e-3;ae5=3.3e-3;ae6=3.3e-3;
ah=3.6e-3;
te=12;
th=100;
Dt=5;
%variable Declaration
Y=zeros(sDIM,6);
X1=zeros(sDIM,1);
%control inputs
filename = 'controller.txt';
input=importdata(filename);
for i=1:length(input)
u(1,i)=bitand(floor(input(i)/2),1);
u(2,i)=bitand(input(i),1);
end
% initial conditions
X_init=[16 16 16 16 16 16];
%X_init=19.5*ones(sDIM,D_tau);
%initial conditions
X=X_init;
time = 0;
k=0;
for j = 1:length(input)*Dt
k=ceil(j/5);
X1(1)=X(1)+Dt*((-3*a-ae1-ah*u(1,k))*X(1)+a*X(2)+a*X(3)+a*X(5)+ae1*te+ah*th*u(1,k));
X1(2)=X(2)+Dt*((-2*a-ae2)*X(2)+a*X(1)+a*X(4)+ae2*te);
X1(3)=X(3)+Dt*((-2*a-ae3)*X(3)+a*X(1)+a*X(4)+ae3*te);
X1(4)=X(4)+Dt*((-3*a-ae4-ah*u(2,k))*X(4)+a*X(2)+a*X(3)+a*X(6)+ae4*te+ah*th*u(2,k));
X1(5)=X(5)+Dt*((-a-ae5)*X(5)+a*X(1)+ae5*te);
X1(6)=X(6)+Dt*((-a-ae6)*X(6)+a*X(4)+ae6*te);
X=X1;
Y(:,j)=X1;
Time(1,j)=time;
time=time+Dt;
end
time1=0;
for l = 1:length(input)*Dt*100
m=ceil(l/500);
Uem(1,l)=u(1,m);
Uem(2,l)=u(2,m);
time1=time1+Dt/100;
Time1(1,l)=time1;
end
figure(1);
plot(Time,Y);
figure(2);
plot(Y(1,:),Y(2,:));
hold on;
plot([17.5, 17.5 ,21,21,17.5],[17.5,21,21,17.5,17.5]);
figure(3);
plot(Y(3,:),Y(4,:));
hold on;
plot([17.5, 17.5 ,21,21,17.5],[17.5,21,21,17.5,17.5]);
figure(4);
plot(Y(5,:),Y(6,:));
hold on;
plot([17.5, 17.5 ,21,21,17.5],[17.5,21,21,17.5,17.5]);
figure(5);
plot(Time1,Uem(1,:));
figure(6);
plot(Time1,Uem(2,:));
#
# compiler
#
#CC = g++
CC = clang++
CXXFLAGS = -Wall -Wextra -std=c++11 -O3 -g -O0
#
# scots
#
SCOTSROOT = ../../..
SCOTSINC = -I$(SCOTSROOT)/bdd -I$(SCOTSROOT)/utils
#
# cudd
#
CUDDPATH = /opt/local
CUDDINC = -I$(CUDDPATH)/include
CUDDLIBS = -lcudd
CUDDLPATH = -L$(CUDDPATH)/lib
TARGET = trafficmodel
all: $(TARGET)
%.o:%.cc
$(CC) -c $(CXXFLAGS) $(CUDDINC) $(SCOTSINC) $< -o $@
$(TARGET): $(TARGET).o
$(CC) $(CXXFLAGS) -o $(TARGET) $(TARGET).o $(CUDDLPATH) $(CUDDLIBS)
clean:
rm ./$(TARGET) ./$(TARGET).o
#include <array>
#include <iostream>
#include<cmath>
#include "cuddObj.hh"
#include "../../src/TicToc.hh"
#include "../../src/fixedPointMode.hh"
#include "../../src/getAbstraction.hh"
#include "../../src/SymbolicSetSpace.hh"
#include "../../src/abstractionMode.hh"
using namespace std;
const int sDIM = 5; /* System dimension */
const int iDIM = 5; /* Input dimension */
size_t P = 3; /* Cardinality of (quantized) input set */
size_t N = 12; /* Temporal Horizon */
typedef std::array<double,sDIM> state_type;
/* ODE */
auto system_post = [](state_type &x, int p) -> void {
state_type z = x;
/* the ode describing the system */
const double v = 70 ;
const double l = 0.25 ;
const double T = 0.002777;
const double q = 0.25;
double b[iDIM];
if(p==0)
{
b[0] = 6;
b[1] = 0;
b[2] = 8;
b[3] = 0;