4341.sagews

t=var("t");
f(t,y) = (7/10)*y*t-(2/10)*t ;
y = function("y", t) ; #defines a variable t and a function y which depends on t

g=desolve(diff(y,t)-f(t,y) , y, ics=[0,1]); g  

#g is the solution to the ODE dy/dt - y - t = 0 with initial condition t=0, y=1. 

print(g); #prints the solution

P_g=plot(g, 0, 4);P_g.show(figsize=4) #plots the solution between t=0 and t=4.

5/7*e^(7/20*t^2) + 2/7

# Compute the next Picard iteration.
# f(t,y) should be the RHS of the DE
# y(t) should be the function to Picard-iterate
# (t0,y0) is the initial value

def PicardIteration (f, y, 0, 1):
    antideriv = integral(f(t,y), t)
    return (1 + antideriv - antideriv(t=0)).simplify()

# define our first iteration to be constantly equal to the initial value
p0(t) = 1           # initial Picard iteration, identically = 1
p = [p0]            # list of Picard iterations (p for Picard)


while len(p)-2 < 10 :  #Computes up to the iterate p_5.
    # how close is each iterate to the actual value at t=2?
    print(len(p)-1, (p[-1])(4), g(4), var('difference') == float( abs((p[-1])(4) - g(4))))
    p.append(PicardIteration(f, p[-1], 0, 1))

(0, 1, 5/7*e^(28/5) + 2/7, difference == 192.44743387582321) (1, 5, 5/7*e^(28/5) + 2/7, difference == 188.44743387582321) (2, 81/5, 5/7*e^(28/5) + 2/7, difference == 177.24743387582322) (3, 2783/75, 5/7*e^(28/5) + 2/7, difference == 156.34076720915655) (4, 8297/125, 5/7*e^(28/5) + 2/7, difference == 127.07143387582322) (5, 929603/9375, 5/7*e^(28/5) + 2/7, difference == 94.289780542489893) (6, 18246637/140625, 5/7*e^(28/5) + 2/7, difference == 63.693570764712121) (7, 36147851/234375, 5/7*e^(28/5) + 2/7, difference == 39.216602942489885) (8, 602454053/3515625, 5/7*e^(28/5) + 2/7, difference == 22.082725466934335) (9, 28797048449/158203125, 5/7*e^(28/5) + 2/7, difference == 11.421646148810879) (10, 247846278707/1318359375, 5/7*e^(28/5) + 2/7, difference == 5.4514417306617418)

G=plot((p[0], p[1], p[2], p[3], p[4], p[5]), 0, 2, rgbcolor=hue(1));
G+=P_g;
G.show(figsize=5)