PYOMO ODE Solve with Manual Discretization

[Audrey Demafo]

Hi, can anyone help me figure out where I made a mistake in the following code ? I am trying to solve a simple ode (dxdt = 5x -3, with initial conditions x(0)=0.2). However, the solution I am getting is wrong as it does not match the values obtained with the analytical solution which (x(t) = -0.4*exp(5*t) + 0.6). 

m = ConcreteModel()
m.tf = Param(initialize=24)
m.t = RangeSet(0, value(m.tf))
m.x = Var(m.t)
m.ht = Param(initialize=(value(m.tf))/(value(m.tf)))

m.dxdt = Var(m.t)

## Define the ODE
def _ode(m, k):
    return m.dxdt[k]==5*m.x[k] - 3

## Discretization scheme
def _ode_discr(m, k):
    if k==0:
        return Constraint.Skip
    return m.x[k] == m.x[k-1] + m.ht*m.dxdt[k-1]   ## Forward euler discretization

def _initial_cond(m):
    return m.x[0]==0.2  
m.initial_cond = Constraint(rule=_initial_cond)   

m.obj = 1

sol=[]

for t in m.t:
    if 'ode' and 'ode_discr' in m.component_map():
        m.del_component(m.ode)
        m.del_component(m.ode_discr)
    m.ode = Constraint(m.t, rule=_ode)
    m.ode_discr = Constraint(m.t, rule=_ode_discr)
    solver = SolverFactory("ipopt")
    solver.solve(m, tee=True)
    print("at", t, "", "value=", value(m.x[t]))
    sol.append(value(m.x[t]))
print(sol)

Any help on this would be greatly appreciated. 

Best regards,

Audrey.

DISCLAIMER: The contents of this email and any attachments are confidential. They are intended for the named recipient(s) only. If you have received this email by mistake, please notify the sender immediately and you are herewith notified that the contents are legally privileged and that you do not have permission to disclose the contents to anyone, make copies thereof, retain or distribute or act upon it by any means, electronically, digitally or in print. The views expressed in this communication may be of a personal nature and not be representative of AIMS-Rwanda and/or any of its Initiatives.