diff --git a/__pycache__/neuralnetwork.cpython-39.pyc b/__pycache__/neuralnetwork.cpython-39.pyc index 8bebeb3..8bf0a8a 100644 Binary files a/__pycache__/neuralnetwork.cpython-39.pyc and b/__pycache__/neuralnetwork.cpython-39.pyc differ diff --git a/__pycache__/two_area.cpython-39.pyc b/__pycache__/two_area.cpython-39.pyc index 2e9fa18..aac18c3 100644 Binary files a/__pycache__/two_area.cpython-39.pyc and b/__pycache__/two_area.cpython-39.pyc differ diff --git a/neuralnetwork.py b/neuralnetwork.py index d9e920d..9d0e227 100644 --- a/neuralnetwork.py +++ b/neuralnetwork.py @@ -32,6 +32,7 @@ class NeuralNetwork: # Description : Takes in the state vector at a given time step and computes the output vector for the next layer output = 1/(1 + np.exp(self.wh.dot(self.X)) ) + diff --git a/two_area.py b/two_area.py index c1f58ab..92923e3 100644 --- a/two_area.py +++ b/two_area.py @@ -5,11 +5,15 @@ import math class TwoAreaPS: - def __init__(self, Tg, Tp, Tt, Kp, T12, a12, R, T, beta1, beta2, yt_1,yt_2,yt_3): + def __init__(self, Tg, Tp, Tt, Kp, T12, a12, R, T, beta1, beta2, yt_1,yt_2,yt_3, yt_1_,yt_2_,yt_3_): - self.yt_1 = yt_1 - self.yt_2 = yt_2 - self.yt_3 = yt_3 + self.yt_1_a1 = yt_1 + self.yt_2_a1 = yt_2 + self.yt_3_a1 = yt_3 + + self.yt_1_a2 = yt_1_ + self.yt_2_a2 = yt_2_ + self.yt_3_a2 = yt_3_ self.Xprev = np.zeros( (7,1) ) self.Y = np.zeros( (2,1) ) @@ -57,8 +61,9 @@ class TwoAreaPS: [0,0,0,0], [0,1/Tg,0,-Kp/Tp] ]) - self.C = np.array( [ [ self.beta1 , 0 , 0 ,1 , 0, 0 , 0 ] ]) - # [ 0 , 0, 0, 1, self.beta2, 0, 0 ] ] ) + self.C = np.array( [ [ self.beta1 , 0 , 0 ,1 , 0, 0 , 0 ] , + [ 0 , 0, 0, 1, self.beta2, 0, 0 ] , + [ 0 , 0, 0, 1, 0, 0, 0 ] ] ) # Calculating Discrete Coefs @@ -71,14 +76,13 @@ class TwoAreaPS: def Output(self,Ut): - self.yt_1 , self.yt_2 , self.yt_3 = self.Y[0,0], self.yt_1, self.yt_2 + self.yt_1_a1 , self.yt_2_a1 , self.yt_3_a1 = self.Y[0,0], self.yt_1_a1, self.yt_2_a1 + self.yt_1_a2 , self.yt_2_a2 , self.yt_3_a2 = self.Y[1,0], self.yt_1_a2, self.yt_2_a2 self.X = np.dot( self.Ad, self.Xprev ) + np.dot( self.Bd, Ut ) self.Y = np.dot( self.C, self.Xprev) - # print("Chooth :" , self.Y) - self.Xprev = self.X if (math.isnan(self.Y[0,0])): diff --git a/twoarea_with_PIcontroller.py b/twoarea_with_PIcontroller.py new file mode 100644 index 0000000..e6d2d86 --- /dev/null +++ b/twoarea_with_PIcontroller.py @@ -0,0 +1,276 @@ +import matplotlib.pyplot as plt +import numpy as np +from numpy.lib.function_base import append + +import warnings +warnings.filterwarnings("ignore") + + +from two_area import * + +import neuralnetwork + + +def y(yd): + + PL1 = 0 + PL2 = 0 + + # rbf = RBF.RBF( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.95) + + nn1 = neuralnetwork.NeuralNetwork( aw = 0.000065, av = 0.022, au = 0.025 , gamma = 0.98) + nn2 = neuralnetwork.NeuralNetwork( aw = 0.000065, av = 0.022, au = 0.025 , gamma = 0.98) + + # av=0.021,au = 0.025, aw = 0.0003asig = 0.01 gamma = 0.9 + + # ( aw = 0.0003, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00474099 + # ( aw = 0.00008, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00037917 Settling Time ~ 70 + # ( aw = 0.00007, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00033166 Settling Time ~ 60 + # ( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00031556 Settling Time ~ 60 + # ( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.98) ----> Steady State Error : 0.00118022 + + + + + + Tg = 0.08 + Tp = 20 + Tt = 0.3 + Kp = 120 + T12 = 0.545/(2*np.pi) + a12 = -1 + R = 5 + T = dt = 1/80 + beta1 = 0.425 + beta2 = 0.425 + + yt_1 = 0 + yt_2 = 0 + yt_3 = 0 + + yt_1_ = 0 + yt_2_ = 0 + yt_3_ = 0 + + System = TwoAreaPS( Tg, Tp, Tt, Kp, T12, a12, R, T, beta1, beta2, yt_1,yt_2,yt_3 , yt_1_ ,yt_2_,yt_3_ ) + + + + initial_states = [ yt_1, yt_2 , yt_3 , yt_1_ , yt_2_ , yt_3_] + + plot_data = {"ut1":[] , "pl1" : [] , "delF1":[] , "ut2":[] , "pl2" : [] , "delF2":[] , "time" : []} + + # KP1 = 1.32579169e-07 + # KI1 = 3.09193550e-04 + # KD1 = 1.91587817e-08 + + # KP2 = 1.18316184e-07 + # KI2 = 2.70082882e-04 + # KD2 = 1.85194296e-08 + + KP1 = -1.35204764e-06 + KI1 = 4.30298026e-04 + KD1 = -5.84913513e-07 + + KP2 = -1.16255511e-06 + KI2 = 2.09594022e-04 + KD2 = -5.96951511e-07 + + + ut_1 = 0 + ut_2 = 0 + + t = 200 + + y=[] + x=[] + + + for i in range(0, int(t/dt) ): + + # print(System.yt_1) + e_t_a1 = 0 - System.yt_1_a1 + del_y_a1 = System.yt_1_a1 - System.yt_2_a1 + del2_y_a1 = System.yt_1_a1 - 2*System.yt_2_a1 + System.yt_3_a1 + + e_t_a2 = 0 - System.yt_1_a2 + del_y_a2 = System.yt_1_a2 - System.yt_2_a2 + del2_y_a2 = System.yt_1_a2 - 2*System.yt_2_a2 + System.yt_3_a2 + + + # ut_1 = ut_1 + 0.00043*e_t - 0.01*del_y - 0*del2_y + ut_1 = ut_1 + KI1*e_t_a1 - KP1*del_y_a1 - KD1*del2_y_a1 + + ut_2 = ut_2 + KI2*e_t_a2 - KP2*del_y_a2 - KD2*del2_y_a2 + + + plot_data["ut1"].append(ut_1) + plot_data["ut2"].append(ut_2) + + # Load Graph 1 + + # if( i*dt >= 40 ): + # PL = 0.2 + # elif ( i*dt >= 10 ): + # PL = 0.8 + # else: + # PL = 0 + + # Load Graph 2 + + # if( i*dt <=10 ): + # PL = 0.2 + # elif ( i*dt <= 20 ): + # PL = 0.3 + # elif ( i*dt <= 30 ): + # PL = 0.4 + # elif ( i*dt <= 40 ): + # PL = 0.5 + # elif ( i*dt <= 50 ): + # PL = 0.6 + # else: + # PL = 0.6 + + # Load Graph 3 + + # if( i*dt <=10 ): + # PL = 0.2 + # elif ( i*dt <= 20 ): + # PL = 0.3 + # elif ( i*dt <= 30 ): + # PL = 0.4 + # elif ( i*dt <= 40 ): + # PL = 0.5 + # elif ( i*dt <= 50 ): + # PL = 0.6 + # elif ( i*dt <= 60 ): + # PL = 0.6 + # elif ( i*dt <= 70 ): + # PL = 0.5 + # elif ( i*dt <= 80 ): + # PL = 0.4 + # elif ( i*dt <= 90 ): + # PL = 0.3 + # elif ( i*dt <= 100 ): + # PL = 0.2 + # else: + # PL = 0.1 + + # Load Graph 4 + + # PL = np.random.normal(0.5,0.1) + + # Load Graph 5 + + if( (i*dt)%10==0 and (i*dt)!=0 ): + PL1 = np.random.normal(0.5,0.1) + PL2 = np.random.normal(0.5,0.1) + + # Load Graph 6 + + # if ( i*dt > 50 and i*dt <130): + # PL1 = 0.5 + + # else: + # PL1 = 0 + + # if ( i*dt > 20 and i*dt <70): + # PL2 = 0.5 + + # else: + # PL2 = 0 + + + plot_data["pl1"].append(PL1) + plot_data["pl2"].append(PL2) + + Ut = [ [ut_1] , [ut_2] , [ PL1] , [PL2] ] + + System.Output(Ut) + + + plot_data["delF1"].append(System.Y[0,0]) + + plot_data["delF2"].append(System.Y[1,0]) + + plot_data["time"].append(i*dt) + + + print("The steady state error of area 1 system is : ", System.Y[0][0]) + + + print("The steady state error of area 2 system is : ", System.Y[1][0]) + + return plot_data,initial_states + + +if __name__=="__main__": + + + yd = [0 for i in range(200*80) ] + + + ## Generate Reference array here + + plot_data,i = y(yd) + + + plt.subplot(2,3,1) + plt.plot(plot_data["time"],plot_data["pl1"], label="Reference Signal") + plt.title( "Load vs Time") + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + + plt.subplot(2,3,2) + plt.plot(plot_data["time"],plot_data["ut1"], label="Reference Signal") + plt.title( "Control Signal vs Time") + plt.ylabel("Control Signal") + plt.xlabel("Time (s)") + + plt.subplot(2,3,3) + plt.plot(plot_data["time"],yd, label="Reference Signal") + plt.plot(plot_data["time"],plot_data["delF1"],label ="Output") + + + plt.title( " Initial States y(t-1) , y(t-2) and y(t-3) are " + str(i[0]) + ", " + str(i[1]) +" and "+ str(i[2]) ) + + plt.legend() + + plt.subplot(2,3,4) + + plt.plot(plot_data["time"],plot_data["pl2"], label="Reference Signal") + plt.title( "Load vs Time") + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + + + + plt.subplot(2,3,5) + plt.plot(plot_data["time"],plot_data["ut2"], label="Reference Signal") + plt.title( "Control Signal vs Time") + plt.ylabel("Control Signal") + plt.xlabel("Time (s)") + + + + + plt.subplot(2,3,6) + plt.plot(plot_data["time"],yd, label="Reference Signal") + plt.plot(plot_data["time"],plot_data["delF2"],label ="Output") + + + + + + + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + plt.show() + + diff --git a/twoarea_with_adaptive_pid.py b/twoarea_with_adaptive_pid.py index ee8e6e8..585608b 100644 --- a/twoarea_with_adaptive_pid.py +++ b/twoarea_with_adaptive_pid.py @@ -2,14 +2,34 @@ import matplotlib.pyplot as plt import numpy as np from numpy.lib.function_base import append +import warnings +warnings.filterwarnings("ignore") + + from two_area import * import RBF - def y(yd): - rbf = RBF.RBF( aw = 0.0003, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) + PL1 = 0 + PL2 = 0 + + # rbf = RBF.RBF( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.95) + + nn1 = RBF.RBF( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.95) + nn2 = RBF.RBF( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.95) + + # av=0.021,au = 0.025, aw = 0.0003asig = 0.01 gamma = 0.9 + + # ( aw = 0.0003, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00474099 + # ( aw = 0.00008, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00037917 Settling Time ~ 70 + # ( aw = 0.00007, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00033166 Settling Time ~ 60 + # ( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00031556 Settling Time ~ 60 + # ( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.98) ----> Steady State Error : 0.00118022 + + + Tg = 0.08 @@ -18,8 +38,8 @@ def y(yd): Kp = 120 T12 = 0.545/(2*np.pi) a12 = -1 - R = 2.4 - T = dt = 1/400 + R = 5 + T = dt = 1/80 beta1 = 0.425 beta2 = 0.425 @@ -27,22 +47,38 @@ def y(yd): yt_2 = 0 yt_3 = 0 - System = TwoAreaPS( Tg, Tp, Tt, Kp, T12, a12, R, T, beta1, beta2, yt_1,yt_2,yt_3 ) + yt_1_ = 0 + yt_2_ = 0 + yt_3_ = 0 + + System = TwoAreaPS( Tg, Tp, Tt, Kp, T12, a12, R, T, beta1, beta2, yt_1,yt_2,yt_3 , yt_1_ ,yt_2_,yt_3_ ) - initial_states = [ yt_1, yt_2 , yt_3] + initial_states = [ yt_1, yt_2 , yt_3 , yt_1_ , yt_2_ , yt_3_] - plot_data = {"ut":[] , "pl" : [] , "delF":[] , 'KI' : [], 'KP' : [] , 'KD' : [] , "time" : []} + plot_data = {"ut1":[] , "pl1" : [] , "delF1":[] , 'KI1' : [], 'KP1' : [] , 'KD1' : [] , + "ut2":[] , "pl2" : [] , "delF2":[] , 'KI2' : [], 'KP2' : [] , 'KD2' : [] , "time" : []} Ki = 0 Kd = 0 Kp = 0 - ut_1 = 0 - t = 100 + nn1.K[0] = 13 + nn1.K[1] = 1132 + nn1.K[2] = 32 + + nn2.K[0] = 13 + nn2.K[1] = 1132 + nn2.K[2] = 32 + + + ut_1 = 0 + ut_2 = 0 + + t = 200 y=[] x=[] @@ -51,38 +87,141 @@ def y(yd): for i in range(0, int(t/dt) ): # print(System.yt_1) - e_t = 0 - System.yt_1 - del_y = System.yt_1 - System.yt_2 - del2_y = System.yt_1 - 2*System.yt_2 + System.yt_3 + e_t_a1 = 0 - System.yt_1_a1 + del_y_a1 = System.yt_1_a1 - System.yt_2_a1 + del2_y_a1 = System.yt_1_a1 - 2*System.yt_2_a1 + System.yt_3_a1 - rbf.X[:,0] = [ e_t , -del_y , -del2_y] - rbf.HiddenLayer() - rbf.OutputLayer() + nn1.X[:,0] = [ e_t_a1 , -del_y_a1 , -del2_y_a1] + nn1.HiddenLayer() + nn1.OutputLayer() + + e_t_a2 = 0 - System.yt_1_a2 + del_y_a2 = System.yt_1_a2 - System.yt_2_a2 + del2_y_a2 = System.yt_1_a2 - 2*System.yt_2_a2 + System.yt_3_a2 + + nn2.X[:,0] = [ e_t_a2 , -del_y_a2 , -del2_y_a2] + nn2.HiddenLayer() + nn2.OutputLayer() # ut_1 = ut_1 + 0.00043*e_t - 0.01*del_y - 0*del2_y - ut_1 = ut_1 + rbf.K[1]*e_t - rbf.K[0]*del_y - rbf.K[2]*del2_y - - plot_data["ut"].append(ut_1) + ut_1 = ut_1 + nn1.K[1]*e_t_a1 - nn1.K[0]*del_y_a1 - nn1.K[2]*del2_y_a1 + ut_2 = ut_2 + nn2.K[1]*e_t_a2 - nn2.K[0]*del_y_a2 - nn2.K[2]*del2_y_a2 + + + plot_data["ut1"].append(ut_1) + plot_data["ut2"].append(ut_2) + # Load Graph 1 + + # if( i*dt >= 40 ): + # PL = 0.2 + # elif ( i*dt >= 10 ): + # PL = 0.8 + # else: + # PL = 0 + + # Load Graph 2 + + # if( i*dt <=10 ): + # PL = 0.2 + # elif ( i*dt <= 20 ): + # PL = 0.3 + # elif ( i*dt <= 30 ): + # PL = 0.4 + # elif ( i*dt <= 40 ): + # PL = 0.5 + # elif ( i*dt <= 50 ): + # PL = 0.6 + # else: + # PL = 0.6 + + # Load Graph 3 + + # if( i*dt <=10 ): + # PL = 0.2 + # elif ( i*dt <= 20 ): + # PL = 0.3 + # elif ( i*dt <= 30 ): + # PL = 0.4 + # elif ( i*dt <= 40 ): + # PL = 0.5 + # elif ( i*dt <= 50 ): + # PL = 0.6 + # elif ( i*dt <= 60 ): + # PL = 0.6 + # elif ( i*dt <= 70 ): + # PL = 0.5 + # elif ( i*dt <= 80 ): + # PL = 0.4 + # elif ( i*dt <= 90 ): + # PL = 0.3 + # elif ( i*dt <= 100 ): + # PL = 0.2 + # else: + # PL = 0.1 + + # Load Graph 4 + + # PL = np.random.normal(0.5,0.1) + + # Load Graph 5 + + if( (i*dt)%10==0 and (i*dt)!=0 ): + PL1 = np.random.normal(0.5,0.1) + PL2 = np.random.normal(0.5,0.1) + + # Load Graph 6 + + # if ( i*dt > 50 and i*dt <130): + # PL1 = 0.5 + + # else: + # PL1 = 0 + + # if ( i*dt > 20 and i*dt <70): + # PL2 = 0.5 + + # else: + # PL2 = 0 - PL = 0.2 if( i*dt >= 0.2 ) else 0 - plot_data["pl"].append(PL) + + plot_data["pl1"].append(PL1) + plot_data["pl2"].append(PL2) - Ut = [ [ut_1] , [0] , [ PL] , [0] ] + Ut = [ [ut_1] , [ut_2] , [ PL1] , [PL2] ] System.Output(Ut) - print(rbf.K) + # print(rbf.K) - rbf.Update(0 ,System.Y[0,0] ,System.yt_1 , System.yt_2, System.yt_3 ) + nn1.Update(0 ,System.Y[0,0] ,System.yt_1_a1 , System.yt_2_a1, System.yt_3_a1 ) + nn2.Update(0 ,System.Y[1,0] ,System.yt_1_a2 , System.yt_2_a2, System.yt_3_a2 ) + + + plot_data["delF1"].append(System.Y[0,0]) + + plot_data["delF2"].append(System.Y[1,0]) + + plot_data["KI1"].append(nn1.K[1]) + plot_data["KP1"].append(nn1.K[0]) + plot_data["KD1"].append(nn1.K[2]) + + + plot_data["KI2"].append(nn2.K[1]) + plot_data["KP2"].append(nn2.K[0]) + plot_data["KD2"].append(nn2.K[2]) - plot_data["delF"].append(System.Y[0,0]) plot_data["time"].append(i*dt) - plot_data["KI"].append(rbf.K[1]) - plot_data["KP"].append(rbf.K[0]) - plot_data["KD"].append(rbf.K[2]) + + print("Final Tuned Kp , Ki and Kd values for area 1 are :" , nn1.K) + + print("The steady state error of area 1 system is : ", System.Y[0][0]) + + print("Final Tuned Kp , Ki and Kd values for area 2 are :" , nn2.K) + + print("The steady state error of area 2 system is : ", System.Y[1][0]) return plot_data,initial_states @@ -90,7 +229,7 @@ def y(yd): if __name__=="__main__": - yd = [0 for i in range(100*400) ] + yd = [0 for i in range(200*80) ] ## Generate Reference array here @@ -99,7 +238,7 @@ if __name__=="__main__": plt.subplot(2,3,1) - plt.plot(plot_data["time"],plot_data["pl"], label="Reference Signal") + plt.plot(plot_data["time"],plot_data["pl1"], label="Reference Signal") plt.title( "Load vs Time") plt.ylabel(" Output from System ") plt.xlabel("Time (s)") @@ -107,14 +246,14 @@ if __name__=="__main__": plt.subplot(2,3,2) - plt.plot(plot_data["time"],plot_data["ut"], label="Reference Signal") + plt.plot(plot_data["time"],plot_data["ut1"], label="Reference Signal") plt.title( "Control Signal vs Time") plt.ylabel("Control Signal") plt.xlabel("Time (s)") plt.subplot(2,3,3) plt.plot(plot_data["time"],yd, label="Reference Signal") - plt.plot(plot_data["time"],plot_data["delF"],label ="Output") + plt.plot(plot_data["time"],plot_data["delF1"],label ="Output") plt.title( " Initial States y(t-1) , y(t-2) and y(t-3) are " + str(i[0]) + ", " + str(i[1]) +" and "+ str(i[2]) ) @@ -122,20 +261,81 @@ if __name__=="__main__": plt.legend() plt.subplot(2,3,4) - plt.plot(plot_data["time"],plot_data["KI"], label="KI") + + plt.plot(plot_data["time"],plot_data["pl2"], label="Reference Signal") + plt.title( "Load vs Time") + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + + + + plt.subplot(2,3,5) + plt.plot(plot_data["time"],plot_data["ut2"], label="Reference Signal") + plt.title( "Control Signal vs Time") + plt.ylabel("Control Signal") + plt.xlabel("Time (s)") + + + + + plt.subplot(2,3,6) + plt.plot(plot_data["time"],yd, label="Reference Signal") + plt.plot(plot_data["time"],plot_data["delF2"],label ="Output") + + + + + + + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + plt.show() + + + + plt.subplot(2,3,1) + plt.plot(plot_data["time"],plot_data["KI1"], label="KI") + plt.title( "KI vs Time") + plt.ylabel("KI") + plt.xlabel("Time (s)") + + plt.subplot(2,3,2) + plt.plot(plot_data["time"],plot_data["KP1"], label="KP") + plt.title( "KP vs Time") + plt.ylabel("KP ") + plt.xlabel("Time (s)") + + plt.subplot(2,3,3) + plt.plot(plot_data["time"],plot_data["KD1"], label="KD") + plt.title( "KD vs Time") + plt.ylabel("KD") + plt.xlabel("Time (s)") + + + + plt.title( " Initial States y(t-1) , y(t-2) and y(t-3) are " + str(i[0]) + ", " + str(i[1]) +" and "+ str(i[2]) ) + + plt.legend() + + plt.subplot(2,3,4) + plt.plot(plot_data["time"],plot_data["KI2"], label="KI") plt.title( "KI vs Time") plt.ylabel("KI") plt.xlabel("Time (s)") plt.subplot(2,3,5) - plt.plot(plot_data["time"],plot_data["KP"], label="KP") + plt.plot(plot_data["time"],plot_data["KP2"], label="KP") plt.title( "KP vs Time") plt.ylabel("KP ") plt.xlabel("Time (s)") plt.subplot(2,3,6) - plt.plot(plot_data["time"],plot_data["KD"], label="KD") + plt.plot(plot_data["time"],plot_data["KD2"], label="KD") plt.title( "KD vs Time") plt.ylabel("KD") plt.xlabel("Time (s)") @@ -148,4 +348,3 @@ if __name__=="__main__": plt.show() - diff --git a/twoarea_with_adaptive_pid_nn.py b/twoarea_with_adaptive_pid_nn.py new file mode 100644 index 0000000..f79badc --- /dev/null +++ b/twoarea_with_adaptive_pid_nn.py @@ -0,0 +1,361 @@ +import matplotlib.pyplot as plt +import numpy as np +from numpy.lib.function_base import append + +import warnings +warnings.filterwarnings("ignore") + + +from two_area import * + +import neuralnetwork + + +def y(yd): + + PL1 = 0 + PL2 = 0 + + # rbf = RBF.RBF( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.95) + + nn1 = neuralnetwork.NeuralNetwork( aw = 0.000065, av = 0.022, au = 0.025 , gamma = 0.98) + nn2 = neuralnetwork.NeuralNetwork( aw = 0.000065, av = 0.022, au = 0.025 , gamma = 0.98) + + # av=0.021,au = 0.025, aw = 0.0003asig = 0.01 gamma = 0.9 + + # ( aw = 0.0003, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00474099 + # ( aw = 0.00008, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00037917 Settling Time ~ 70 + # ( aw = 0.00007, av = 0.021, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00033166 Settling Time ~ 60 + # ( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.9) ---> Steady State Error : 0.00031556 Settling Time ~ 60 + # ( aw = 0.000065, av = 0.022, au = 0.025 , asig = 0.01, gamma = 0.98) ----> Steady State Error : 0.00118022 + + + + + + Tg = 0.08 + Tp = 20 + Tt = 0.3 + Kp = 120 + T12 = 0.545/(2*np.pi) + a12 = -1 + R = 5 + T = dt = 1/80 + beta1 = 0.425 + beta2 = 0.425 + + yt_1 = 0 + yt_2 = 0 + yt_3 = 0 + + yt_1_ = 0 + yt_2_ = 0 + yt_3_ = 0 + + System = TwoAreaPS( Tg, Tp, Tt, Kp, T12, a12, R, T, beta1, beta2, yt_1,yt_2,yt_3 , yt_1_ ,yt_2_,yt_3_ ) + + + + initial_states = [ yt_1, yt_2 , yt_3 , yt_1_ , yt_2_ , yt_3_] + + plot_data = {"ut1":[] , "pl1" : [] , "delF1":[] , 'KI1' : [], 'KP1' : [] , 'KD1' : [] , + "ut2":[] , "pl2" : [] , "delF2":[] , 'KI2' : [], 'KP2' : [] , 'KD2' : [] , 'Tie Line' : [ ] ,"time" : []} + + + Ki = 0 + Kd = 0 + Kp = 0 + + + nn1.K[0] = 13 + nn1.K[1] = 1132 + nn1.K[2] = 32 + + nn2.K[0] = 13 + nn2.K[1] = 1132 + nn2.K[2] = 32 + + + ut_1 = 0 + ut_2 = 0 + + t = 200 + + y=[] + x=[] + + + for i in range(0, int(t/dt) ): + + # print(System.yt_1) + e_t_a1 = 0 - System.yt_1_a1 + del_y_a1 = System.yt_1_a1 - System.yt_2_a1 + del2_y_a1 = System.yt_1_a1 - 2*System.yt_2_a1 + System.yt_3_a1 + + nn1.X[:,0] = [ e_t_a1 , -del_y_a1 , -del2_y_a1] + nn1.HiddenLayer() + nn1.OutputLayer() + + e_t_a2 = 0 - System.yt_1_a2 + del_y_a2 = System.yt_1_a2 - System.yt_2_a2 + del2_y_a2 = System.yt_1_a2 - 2*System.yt_2_a2 + System.yt_3_a2 + + nn2.X[:,0] = [ e_t_a2 , -del_y_a2 , -del2_y_a2] + nn2.HiddenLayer() + nn2.OutputLayer() + + + # ut_1 = ut_1 + 0.00043*e_t - 0.01*del_y - 0*del2_y + ut_1 = ut_1 + nn1.K[1]*e_t_a1 - nn1.K[0]*del_y_a1 - nn1.K[2]*del2_y_a1 + + ut_2 = ut_2 + nn2.K[1]*e_t_a2 - nn2.K[0]*del_y_a2 - nn2.K[2]*del2_y_a2 + + + plot_data["ut1"].append(ut_1) + plot_data["ut2"].append(ut_2) + + # Load Graph 1 + + # if( i*dt >= 40 ): + # PL = 0.2 + # elif ( i*dt >= 10 ): + # PL = 0.8 + # else: + # PL = 0 + + # Load Graph 2 + + # if( i*dt <=10 ): + # PL = 0.2 + # elif ( i*dt <= 20 ): + # PL = 0.3 + # elif ( i*dt <= 30 ): + # PL = 0.4 + # elif ( i*dt <= 40 ): + # PL = 0.5 + # elif ( i*dt <= 50 ): + # PL = 0.6 + # else: + # PL = 0.6 + + # Load Graph 3 + + # if( i*dt <=10 ): + # PL = 0.2 + # elif ( i*dt <= 20 ): + # PL = 0.3 + # elif ( i*dt <= 30 ): + # PL = 0.4 + # elif ( i*dt <= 40 ): + # PL = 0.5 + # elif ( i*dt <= 50 ): + # PL = 0.6 + # elif ( i*dt <= 60 ): + # PL = 0.6 + # elif ( i*dt <= 70 ): + # PL = 0.5 + # elif ( i*dt <= 80 ): + # PL = 0.4 + # elif ( i*dt <= 90 ): + # PL = 0.3 + # elif ( i*dt <= 100 ): + # PL = 0.2 + # else: + # PL = 0.1 + + # Load Graph 4 + + # PL = np.random.normal(0.5,0.1) + + # Load Graph 5 + + # if( (i*dt)%10==0 and (i*dt)!=0 ): + # PL1 = np.random.normal(0.5,0.1) + # PL2 = np.random.normal(0.5,0.1) + + # Load Graph 6 + + if ( i*dt > 50 and i*dt <130): + PL1 = 0.5 + + else: + PL1 = 0 + + if ( i*dt > 20 and i*dt <70): + PL2 = 0.5 + + else: + PL2 = 0 + + + plot_data["pl1"].append(PL1) + plot_data["pl2"].append(PL2) + + Ut = [ [ut_1] , [ut_2] , [ PL1] , [PL2] ] + + System.Output(Ut) + + # print(rbf.K) + + nn1.Update(0 ,System.Y[0,0] ,System.yt_1_a1 , System.yt_2_a1, System.yt_3_a1 ) + nn2.Update(0 ,System.Y[1,0] ,System.yt_1_a2 , System.yt_2_a2, System.yt_3_a2 ) + + + plot_data["delF1"].append(System.Y[0,0]) + + plot_data["delF2"].append(System.Y[1,0]) + + plot_data["KI1"].append(nn1.K[1]) + plot_data["KP1"].append(nn1.K[0]) + plot_data["KD1"].append(nn1.K[2]) + + + plot_data["KI2"].append(nn2.K[1]) + plot_data["KP2"].append(nn2.K[0]) + plot_data["KD2"].append(nn2.K[2]) + + plot_data["time"].append(i*dt) + + plot_data["Tie Line"].append(System.Y[2,0]) + + print("Final Tuned Kp , Ki and Kd values for area 1 are :" , nn1.K) + + print("The steady state error of area 1 system is : ", System.Y[0][0]) + + print("Final Tuned Kp , Ki and Kd values for area 2 are :" , nn2.K) + + print("The steady state error of area 2 system is : ", System.Y[1][0]) + + return plot_data,initial_states + + +if __name__=="__main__": + + + yd = [0 for i in range(200*80) ] + + + ## Generate Reference array here + + plot_data,i = y(yd) + + + plt.subplot(2,3,1) + plt.plot(plot_data["time"],plot_data["pl1"], label="Reference Signal") + plt.title( "Load vs Time") + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + + plt.subplot(2,3,2) + plt.plot(plot_data["time"],plot_data["ut1"], label="Reference Signal") + plt.title( "Control Signal vs Time") + plt.ylabel("Control Signal") + plt.xlabel("Time (s)") + + plt.subplot(2,3,3) + plt.plot(plot_data["time"],yd, label="Reference Signal") + plt.plot(plot_data["time"],plot_data["delF1"],label ="Output") + + + plt.title( " Initial States y(t-1) , y(t-2) and y(t-3) are " + str(i[0]) + ", " + str(i[1]) +" and "+ str(i[2]) ) + + plt.legend() + + plt.subplot(2,3,4) + + plt.plot(plot_data["time"],plot_data["pl2"], label="Reference Signal") + plt.title( "Load vs Time") + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + + + + plt.subplot(2,3,5) + plt.plot(plot_data["time"],plot_data["ut2"], label="Reference Signal") + plt.title( "Control Signal vs Time") + plt.ylabel("Control Signal") + plt.xlabel("Time (s)") + + + + + plt.subplot(2,3,6) + plt.plot(plot_data["time"],yd, label="Reference Signal") + plt.plot(plot_data["time"],plot_data["delF2"],label ="Output") + + + + + + + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + plt.show() + + + + plt.subplot(2,3,1) + plt.plot(plot_data["time"],plot_data["KI1"], label="KI") + plt.title( "KI vs Time") + plt.ylabel("KI") + plt.xlabel("Time (s)") + + plt.subplot(2,3,2) + plt.plot(plot_data["time"],plot_data["KP1"], label="KP") + plt.title( "KP vs Time") + plt.ylabel("KP ") + plt.xlabel("Time (s)") + + plt.subplot(2,3,3) + plt.plot(plot_data["time"],plot_data["KD1"], label="KD") + plt.title( "KD vs Time") + plt.ylabel("KD") + plt.xlabel("Time (s)") + + + + plt.title( " Initial States y(t-1) , y(t-2) and y(t-3) are " + str(i[0]) + ", " + str(i[1]) +" and "+ str(i[2]) ) + + plt.legend() + + plt.subplot(2,3,4) + plt.plot(plot_data["time"],plot_data["KI2"], label="KI") + plt.title( "KI vs Time") + plt.ylabel("KI") + plt.xlabel("Time (s)") + + + plt.subplot(2,3,5) + plt.plot(plot_data["time"],plot_data["KP2"], label="KP") + plt.title( "KP vs Time") + plt.ylabel("KP ") + plt.xlabel("Time (s)") + + plt.subplot(2,3,6) + plt.plot(plot_data["time"],plot_data["KD2"], label="KD") + plt.title( "KD vs Time") + plt.ylabel("KD") + plt.xlabel("Time (s)") + + + + + plt.ylabel(" Output from System ") + plt.xlabel("Time (s)") + + + plt.show() + + plt.plot(plot_data["time"],plot_data["Tie Line"]) + plt.title( "Tie Line vs Time") + plt.ylabel("Tie Line") + plt.xlabel("Time (s)") + + + plt.show()