forwardDQSS_expression.sagews

#Evaluate all cells via pressing the Shift+Return keys or the small "play" button in the menu.
#Newton polynomial: https://en.wikipedia.org/wiki/Newton_polynomial

These are the expression obtained when doing direct replacement in the equations

print("---- Forward QSS1 ----")

# qss1
from sage.symbolic.function_factory import function

var('t y0 ty r0 tr d0  td  e')

def evTau(self, t): return r0
tau = function("tau", nargs=1, eval_func=evTau)
def evY(self, t): return y0
y = function("y", nargs=1, eval_func=evY)
def evD(self, t): return d0
d = function("d", nargs=1, eval_func=evD)

# forward
forward = d(t + tau(t))
expression = (d(t + tau(t)) ==  y(t))


print("---- Forward d(t + tau(t)) ----")
print("Forward d(t + tau(t)) = %s" % forward.simplify())
print("Forward d(t + tau(t)) = %s" % forward.simplify_full())
print("expression = %s" % expression.simplify_full())
print("---- coefficients:----")
coeff0 = expression.substitute(t = 0) # .expand().collect(t)
print("do = %s" % solve(coeff0, d0)[0].right_hand_side().simplify_full()) #  .right_hand_side().simplify_full()

---- Forward QSS1 ---- (t, y0, ty, r0, tr, d0, td, e) ---- Forward d(t + tau(t)) ---- Forward d(t + tau(t)) = d0 Forward d(t + tau(t)) = d0 expression = d0 == y0 ---- coefficients:---- do = y0

print("---- Forward QSS2 ----")


# qss2
from sage.symbolic.function_factory import function

var('t y0 y1 ty r0 r1 tr d0 d1 td  e')

def evTau(self, t): return r0 + r1* (t - tr)
tau = function("tau", nargs=1, eval_func=evTau)
def evY(self, t): return y0 + y1* (t - ty)
y = function("y", nargs=1, eval_func=evY)
def evD(self, t): return d0 + d1* (t - td)
d = function("d", nargs=1, eval_func=evD)

# forward
forward = d(t + tau(t))
expression = (d(t + tau(t)) ==  y(t))


print("---- Forward d(t + tau(t)) ----")
print("Forward d(t + tau(t)) = %s" % forward.simplify())
print("Forward d(t + tau(t)) = %s" % forward.simplify_full())
print("expression = %s" % expression.simplify_full())
print("---- coefficients:----")
coeff0 = expression.substitute(t = 0) # .expand().collect(t)
coeff0 = solve(coeff0, d0)[0].right_hand_side()
print("do = %s \n" % coeff0) #  .right_hand_side().simplify_full()

coeff1 = expression.substitute(d0 = coeff0) # .expand().collect(t)
coeff1 = solve(coeff1, d1)[0].right_hand_side()
print("d1 = %s \n" % coeff1 )

print("---- When new tau arrives:----")
coeff0_newTau = coeff0.substitute(td = tr+r0, d1=coeff1)
print("do = %s (same as the definition of y(tr)=%s)\n" %  (coeff0_newTau.simplify_full(), y(tr) )) #

print("---- When new y arrives:----")
coeff0_newTau = coeff0.substitute(td = ty+tau(ty), d1=coeff1)
print("do = %s (same as the definition of y(ty)=%s)\n" %  (coeff0_newTau.simplify_full(), y(ty) )) #

# print("---- coefficients: in case QSS1 (y1=0, r1=0)----")
# expression = expression.substitute(y1=0, r1=0)
# coeff0 = expression.substitute(t = 0) # .expand().collect(t)
# print("do = %s" % solve(coeff0, d0)[0].right_hand_side().simplify_full()) #  .right_hand_side().simplify_full()
# #coeff1 = expression.substitute(d0 = 0, y0 = 0, r0 = 0)
# coeff1 = expression - coeff0
# print("d1 = %s" % solve(coeff1, d1)[0].right_hand_side().simplify_full())

---- Forward QSS2 ---- (t, y0, y1, ty, r0, r1, tr, d0, d1, td, e) ---- Forward d(t + tau(t)) ---- Forward d(t + tau(t)) = (r1*(t - tr) + r0 + t - td)*d1 + d0 Forward d(t + tau(t)) = -d1*r1*tr + d1*r0 + (d1*r1 + d1)*t - d1*td + d0 expression = -d1*r1*tr + d1*r0 + (d1*r1 + d1)*t - d1*td + d0 == (t - ty)*y1 + y0 ---- coefficients:---- do = d1*r1*tr - d1*r0 + d1*td - ty*y1 + y0 d1 = y1/(r1 + 1) ---- When new tau arrives:---- do = (tr - ty)*y1 + y0 (same as the definition of y(tr)=(tr - ty)*y1 + y0) ---- When new y arrives:---- do = y0 (same as the definition of y(ty)=y0)

print("---- Forward QSS3 ----")
# qss3
from sage.symbolic.function_factory import function

var('t T_y T_r T_d y_0 y_1 y_2 r_0 r_1 r_2 d_0 d_1 d_2')

def evTau(self, t): return r_0 + r_1* (t - T_r) + r_2* (t - T_r)^2
tau = function("tau", nargs=1, eval_func=evTau)
def evY(self, t): return y_0 + y_1* (t - T_y) + y_2* (t - T_y)^2
y = function("y", nargs=1, eval_func=evY)
def evD(self, t): return d_0 + d_1* (t - T_d) + d_2* (t - T_d)^2
d = function("d", nargs=1, eval_func=evD)


# forward
forward = d(t + tau(t))

print("---- Forward d(t + tau(t)) ----")
print("Forward d(t + tau(t)) = %s \n" % forward.simplify())
print("Forward d(t + tau(t)) simplified = %s \n" % forward.simplify_full())

#forward = d_2*r2^2*t^4 + d_2*r2^2*tr^4 + 2*(d_2*r1 + d_2)*r2*t^3 - 2*(2*d_2*r2^2*t + d_2*r1*r2)*tr^3 + d_2*r0^2 + (d_2*r1^2 + 2*d_2*r1 + (2*d_2*r0 + d1)*r2 + d_2)*t^2 + d_2*td^2 + (6*d_2*r2^2*t^2 + d_2*r1^2 + 2*(3*d_2*r1 + d_2)*r2*t - 2*d_2*r2*td + (2*d_2*r0 + d1)*r2)*tr^2 + d1*r0 + (2*d_2*r0 + (2*d_2*r0 + d1)*r1 + d1)*t - (2*d_2*r2*t^2 + 2*d_2*r0 + 2*(d_2*r1 + d_2)*t + d1)*td - (4*d_2*r2^2*t^3 + 2*(3*d_2*r1 + 2*d_2)*r2*t^2 + (2*d_2*r0 + d1)*r1 + 2*(d_2*r1^2 + d_2*r1 + (2*d_2*r0 + d1)*r2)*t - 2*(2*d_2*r2*t + d_2*r1)*td)*tr + d0


# remove t^4 and t^3 (manually :-()
print("---- Removing 3 and 4 order ----")

forward = T_r^4*d_2*r_2^2 + T_r^2*d_2*r_1^2 + T_d^2*d_2 + d_2*r_0^2 + (6*T_r^2*d_2*r_2^2 + d_2*r_1^2 + 2*d_2*r_1 - (6*T_r*d_2*r_1 + 2*(T_d + 2*T_r)*d_2 - 2*d_2*r_0 - d_1)*r_2 + d_2)*t^2 - T_d*d_1 - (2*T_d*d_2 - d_1)*r_0 + (2*T_d*T_r*d_2 - 2*T_r*d_2*r_0 - T_r*d_1)*r_1 - (2*T_r^3*d_2*r_1 + 2*T_d*T_r^2*d_2 - 2*T_r^2*d_2*r_0 - T_r^2*d_1)*r_2 - (4*T_r^3*d_2*r_2^2 + 2*T_r*d_2*r_1^2 + 2*T_d*d_2 - 2*d_2*r_0 + (2*(T_d + T_r)*d_2 - 2*d_2*r_0 - d_1)*r_1 - 2*(3*T_r^2*d_2*r_1 - 2*T_r*d_2*r_0 - T_r*d_1 + (2*T_d*T_r + T_r^2)*d_2)*r_2 - d_1)*t + d_0

print("d(t + tau(t)) = %s\n" % forward)
print("y(t) = %s \n" % y(t))
#
# results in
# t^0:
# T_y^2*y_2  +  - T_y*y_1  + y_0 = + d_0 + T_r^4*d_2*r_2^2 + T_r^2*d_2*r_1^2  + T_d^2*d_2 + d_2*r_0^2  - T_d*d_1 - (2*T_d*d_2 - d_1)*r_0 + (2*T_d*T_r*d_2 - 2*T_r*d_2*r_0 - T_r*d_1)*r_1 - (2*T_r^3*d_2*r_1 + 2*T_d*T_r^2*d_2 - 2*T_r^2*d_2*r_0 - T_r^2*d_1)*r_2
# t^1:
# t* (y_1 - 2*T_y*y_2) = - t * (4*T_r^3*d_2*r_2^2 + 2*T_r*d_2*r_1^2 + 2*T_d*d_2 - 2*d_2*r_0 + (2*(T_d + T_r)*d_2 - 2*d_2*r_0 - d_1)*r_1 - 2*(3*T_r^2*d_2*r_1 - 2*T_r*d_2*r_0 - T_r*d_1 + (2*T_d*T_r + T_r^2)*d_2)*r_2 - d_1)
# t^2:
#t^2*y_2 = t^2 * (6*T_r^2*d_2*r_2^2 + d_2*r_1^2 + 2*d_2*r_1 - (6*T_r*d_2*r_1 + 2*(T_d + 2*T_r)*d_2 - 2*d_2*r_0 - d_1)*r_2 + d_2)


coeff0_y = T_y^2*y_2  +  - T_y*y_1  + y_0
coeff0_d = d_0 + T_r^4*d_2*r_2^2 + T_r^2*d_2*r_1^2  + T_d^2*d_2 + d_2*r_0^2  - T_d*d_1 - (2*T_d*d_2 - d_1)*r_0 + (2*T_d*T_r*d_2 - 2*T_r*d_2*r_0 - T_r*d_1)*r_1 - (2*T_r^3*d_2*r_1 + 2*T_d*T_r^2*d_2 - 2*T_r^2*d_2*r_0 - T_r^2*d_1)*r_2
expression = (coeff0_y ==  coeff0_d)
print("expression = %s" % expression)
coeff0 = solve(expression, d_0)[0].right_hand_side()
print("do = %s \n" % coeff0.simplify_full()) #  .right_hand_side().simplify_full()


coeff1_y = y_1 - 2*T_y*y_2
coeff1_d = 4*T_r^3*d_2*r_2^2 + 2*T_r*d_2*r_1^2 + 2*T_d*d_2 - 2*d_2*r_0 + (2*(T_d + T_r)*d_2 - 2*d_2*r_0 - d_1)*r_1 - 2*(3*T_r^2*d_2*r_1 - 2*T_r*d_2*r_0 - T_r*d_1 + (2*T_d*T_r + T_r^2)*d_2)*r_2 - d_1
expression = (coeff1_y ==  coeff1_d)
coeff1 = solve(expression, d_1)[0].right_hand_side()
print("d1 = %s \n" % coeff1.simplify_full()) #  .right_hand_side().simplify_full()


coeff2_y = y_2
coeff2_d = 6*T_r^2*d_2*r_2^2 + d_2*r_1^2 + 2*d_2*r_1 - (6*T_r*d_2*r_1 + 2*(T_d + 2*T_r)*d_2 - 2*d_2*r_0 - d_1)*r_2 + d_2
expression = (coeff2_y ==  coeff2_d)
expression = expression.substitute(d_1 = coeff1)
coeff2 = solve(expression, d_2)[0].right_hand_side()
print("d2 = %s \n" % coeff2.simplify_full()) #  .right_hand_side().simplify_full()


print("---- When new tau arrives:----")
coeff0_newTau = coeff0.substitute(T_d = T_r+r_0, d_1=coeff1, d_2=coeff2)
print("do = %s )\n" %  coeff0_newTau.simplify_full()) #
print("y(T_r) = %s \n" %  y(T_r).simplify_full()) #


# expression = (forward ==  y(t))

# print("expression = %s" % expression)
# print("expression = %s" % y(t).expand())

# print("---- coefficients:----")
# coeff0 = expression.substitute(t = 0) # .expand().collect(t)
# coeff0 = solve(coeff0, d_0)[0].right_hand_side()
# print("do = %s \n" % coeff0.simplify_full()) #  .right_hand_side().simplify_full()

# coeff1 = expression.substitute(d_0 = coeff0) # .expand().collect(t)
# coeff1 = solve(coeff1, d_1)[0].right_hand_side()
# print("d1 = %s \n" % coeff1 )

# coeff2 = expression.substitute(d_1 = coeff1) # .expand().collect(t)
# coeff2 = solve(coeff2, d_2)[0].right_hand_side()
# print("d2 = %s \n"  % coeff2.simplify_full() )

---- Forward QSS3 ---- (t, T_y, T_r, T_d, y_0, y_1, y_2, r_0, r_1, r_2, d_0, d_1, d_2) ---- Forward d(t + tau(t)) ---- Forward d(t + tau(t)) = ((T_r - t)^2*r_2 - (T_r - t)*r_1 - T_d + r_0 + t)^2*d_2 + ((T_r - t)^2*r_2 - (T_r - t)*r_1 - T_d + r_0 + t)*d_1 + d_0 Forward d(t + tau(t)) simplified = T_r^4*d_2*r_2^2 + d_2*r_2^2*t^4 + T_r^2*d_2*r_1^2 - 2*(2*T_r*d_2*r_2^2 - (d_2*r_1 + d_2)*r_2)*t^3 + T_d^2*d_2 + d_2*r_0^2 + (6*T_r^2*d_2*r_2^2 + d_2*r_1^2 + 2*d_2*r_1 - (6*T_r*d_2*r_1 + 2*(T_d + 2*T_r)*d_2 - 2*d_2*r_0 - d_1)*r_2 + d_2)*t^2 - T_d*d_1 - (2*T_d*d_2 - d_1)*r_0 + (2*T_d*T_r*d_2 - 2*T_r*d_2*r_0 - T_r*d_1)*r_1 - (2*T_r^3*d_2*r_1 + 2*T_d*T_r^2*d_2 - 2*T_r^2*d_2*r_0 - T_r^2*d_1)*r_2 - (4*T_r^3*d_2*r_2^2 + 2*T_r*d_2*r_1^2 + 2*T_d*d_2 - 2*d_2*r_0 + (2*(T_d + T_r)*d_2 - 2*d_2*r_0 - d_1)*r_1 - 2*(3*T_r^2*d_2*r_1 - 2*T_r*d_2*r_0 - T_r*d_1 + (2*T_d*T_r + T_r^2)*d_2)*r_2 - d_1)*t + d_0 ---- Removing 3 and 4 order ---- d(t + tau(t)) = T_r^4*d_2*r_2^2 + T_r^2*d_2*r_1^2 + T_d^2*d_2 + d_2*r_0^2 + (6*T_r^2*d_2*r_2^2 + d_2*r_1^2 + 2*d_2*r_1 - (6*T_r*d_2*r_1 + 2*(T_d + 2*T_r)*d_2 - 2*d_2*r_0 - d_1)*r_2 + d_2)*t^2 - T_d*d_1 - (2*T_d*d_2 - d_1)*r_0 + (2*T_d*T_r*d_2 - 2*T_r*d_2*r_0 - T_r*d_1)*r_1 - (2*T_r^3*d_2*r_1 + 2*T_d*T_r^2*d_2 - 2*T_r^2*d_2*r_0 - T_r^2*d_1)*r_2 - (4*T_r^3*d_2*r_2^2 + 2*T_r*d_2*r_1^2 + 2*T_d*d_2 - 2*d_2*r_0 + (2*(T_d + T_r)*d_2 - 2*d_2*r_0 - d_1)*r_1 - 2*(3*T_r^2*d_2*r_1 - 2*T_r*d_2*r_0 - T_r*d_1 + (2*T_d*T_r + T_r^2)*d_2)*r_2 - d_1)*t + d_0 y(t) = (T_y - t)^2*y_2 - (T_y - t)*y_1 + y_0 expression = T_y^2*y_2 - T_y*y_1 + y_0 == T_r^4*d_2*r_2^2 + T_r^2*d_2*r_1^2 + T_d^2*d_2 + d_2*r_0^2 - T_d*d_1 - (2*T_d*d_2 - d_1)*r_0 + (2*T_d*T_r*d_2 - 2*T_r*d_2*r_0 - T_r*d_1)*r_1 - (2*T_r^3*d_2*r_1 + 2*T_d*T_r^2*d_2 - 2*T_r^2*d_2*r_0 - T_r^2*d_1)*r_2 + d_0 do = -T_r^4*d_2*r_2^2 - T_r^2*d_2*r_1^2 - T_d^2*d_2 - d_2*r_0^2 + T_y^2*y_2 + T_d*d_1 + (2*T_d*d_2 - d_1)*r_0 - (2*T_d*T_r*d_2 - 2*T_r*d_2*r_0 - T_r*d_1)*r_1 + (2*T_r^3*d_2*r_1 + 2*T_d*T_r^2*d_2 - 2*T_r^2*d_2*r_0 - T_r^2*d_1)*r_2 - T_y*y_1 + y_0 d1 = -(4*T_r^3*d_2*r_2^2 + 2*T_r*d_2*r_1^2 + 2*T_d*d_2 - 2*d_2*r_0 + 2*((T_d + T_r)*d_2 - d_2*r_0)*r_1 - 2*(3*T_r^2*d_2*r_1 - 2*T_r*d_2*r_0 + (2*T_d*T_r + T_r^2)*d_2)*r_2 + 2*T_y*y_2 - y_1)/(2*T_r*r_2 - r_1 - 1) d2 = -(r_2*y_1 - (2*(T_r + T_y)*r_2 - r_1 - 1)*y_2)/(8*T_r^3*r_2^3 - r_1^3 - 12*(T_r^2*r_1 + T_r^2)*r_2^2 - 3*r_1^2 + 6*(T_r*r_1^2 + 2*T_r*r_1 + T_r)*r_2 - 3*r_1 - 1) ---- When new tau arrives:---- do = (16*T_r^7*d_2*r_2^5 - 2*T_r^2*d_2*r_1^5 + 2*(T_r*d_2*r_0 - (T_d*T_r + 4*T_r^2)*d_2)*r_1^4 - 8*(7*T_r^6*d_2*r_1 - 2*T_r^5*d_2*r_0 + (2*T_d*T_r^5 + 5*T_r^6)*d_2)*r_2^4 + 4*(2*T_r*d_2*r_0 - (2*T_d*T_r + 3*T_r^2)*d_2)*r_1^3 + 4*(19*T_r^5*d_2*r_1^2 - 10*T_r^4*d_2*r_0 + (10*T_d*T_r^4 + 9*T_r^5)*d_2 - 2*(5*T_r^4*d_2*r_0 - (5*T_d*T_r^4 + 14*T_r^5)*d_2)*r_1)*r_2^3 - 2*T_d*T_r*d_2 + 2*T_r*d_2*r_0 + 4*(3*T_r*d_2*r_0 - (3*T_d*T_r + 2*T_r^2)*d_2)*r_1^2 - 2*(25*T_r^4*d_2*r_1^3 - 18*T_r^3*d_2*r_0 - 3*(6*T_r^3*d_2*r_0 - (6*T_d*T_r^3 + 19*T_r^4)*d_2)*r_1^2 + (18*T_d*T_r^3 + 7*T_r^4)*d_2 - 3*(12*T_r^3*d_2*r_0 - (12*T_d*T_r^3 + 13*T_r^4)*d_2)*r_1)*r_2^2 + 2*(4*T_r*d_2*r_0 - (4*T_d*T_r + T_r^2)*d_2)*r_1 + 2*(8*T_r^3*d_2*r_1^4 - 7*T_r^2*d_2*r_0 - (7*T_r^2*d_2*r_0 - (7*T_d*T_r^2 + 25*T_r^3)*d_2)*r_1^3 - 3*(7*T_r^2*d_2*r_0 - (7*T_d*T_r^2 + 9*T_r^3)*d_2)*r_1^2 + (7*T_d*T_r^2 + T_r^3)*d_2 - (21*T_r^2*d_2*r_0 - (21*T_d*T_r^2 + 11*T_r^3)*d_2)*r_1)*r_2 + (8*T_r^3*r_2^3 - r_1^3 - 12*(T_r^2*r_1 + T_r^2)*r_2^2 - 3*r_1^2 + 6*(T_r*r_1^2 + 2*T_r*r_1 + T_r)*r_2 - 3*r_1 - 1)*y_0 + ((T_r + T_y)*r_1^3 - (3*T_r^4 + 8*T_r^3*T_y)*r_2^3 + 3*(T_r + T_y)*r_1^2 + 6*(T_r^3 + 2*T_r^2*T_y + (T_r^3 + 2*T_r^2*T_y)*r_1)*r_2^2 + 3*(T_r + T_y)*r_1 - 2*((2*T_r^2 + 3*T_r*T_y)*r_1^2 + 2*T_r^2 + 3*T_r*T_y + 2*(2*T_r^2 + 3*T_r*T_y)*r_1)*r_2 + T_r + T_y)*y_1 + ((T_r^2 - 2*T_r*T_y - T_y^2)*r_1^3 - 2*(T_r^5 - 3*T_r^4*T_y - 4*T_r^3*T_y^2)*r_2^3 + 3*(T_r^2 - 2*T_r*T_y - T_y^2)*r_1^2 + (5*T_r^4 - 12*T_r^3*T_y - 12*T_r^2*T_y^2 + (5*T_r^4 - 12*T_r^3*T_y - 12*T_r^2*T_y^2)*r_1)*r_2^2 + T_r^2 - 2*T_r*T_y - T_y^2 + 3*(T_r^2 - 2*T_r*T_y - T_y^2)*r_1 - 2*(2*T_r^3 - 4*T_r^2*T_y - 3*T_r*T_y^2 + (2*T_r^3 - 4*T_r^2*T_y - 3*T_r*T_y^2)*r_1^2 + 2*(2*T_r^3 - 4*T_r^2*T_y - 3*T_r*T_y^2)*r_1)*r_2)*y_2)/(8*T_r^3*r_2^3 - r_1^3 - 12*(T_r^2*r_1 + T_r^2)*r_2^2 - 3*r_1^2 + 6*(T_r*r_1^2 + 2*T_r*r_1 + T_r)*r_2 - 3*r_1 - 1) ) y(T_r) = (T_r - T_y)*y_1 + (T_r^2 - 2*T_r*T_y + T_y^2)*y_2 + y_0

print("---- test:----")
from sage.symbolic.expression_conversions import PolynomialConverter
# to poly
c1 = forward
m=c1.list()
poly = c1.polynomial(SR)

print("c1 = %s" % c1)
print("m1 = %s" % m)
print("poly = %s" % poly)
vars = poly.variables()
poly.dict()

poly[1]
z = SR[vars]
z
# z.dict()
# polynomial(vars, coeffs)

print("---- coefficients: in case QSS2 (y2=0, r2=0)----")
expression = expression.substitute(y2=0, r2=0, d_2=0)
coeff0 = expression.substitute(t = 0) # .expand().collect(t)
coeff0 = solve(coeff0, d_0)[0].right_hand_side()
print("do = %s \n" % coeff0.simplify_full()) #  .right_hand_side().simplify_full()

coeff1 = expression.substitute(d_0 = coeff0) # .expand().collect(t)
coeff1 = solve(coeff1, d_1)[0].right_hand_side()
print("d1 = %s \n" % coeff1.simplify_full())


# print("---- coefficients: in case QSS1 (y2=0, r2=0, y1=0, r1=0)----")
# expression = expression.substitute(y1=0, r1=0, y2=0, r2=0)
# coeff0 = expression.substitute(t = 0) # .expand().collect(t)
# print("do = %s" % solve(coeff0, d0)[0].right_hand_side().simplify_full()) #  .right_hand_side().simplify_full()
# #coeff1 = expression.substitute(d0 = 0, y0 = 0, r0 = 0)
# coeff1 = expression - coeff0
# print("d1 = %s" % solve(coeff1, d1)[0].right_hand_side().simplify_full())