241.sagews

#
# Euler's method #2
#
# dy/dx = 0.2 * x * y
#

from pylab import *

# time step
h = 0.1

# initial condition
x0 = 1.0
y0 = 1.0

# final point
xn = 5.5

# function (DE)
def fn(x, y):
	return 0.2 * x * y;

# exact solution
def fn_exact(x):
	return exp(0.1*(x*x - 1))

# algorithm
y = y0
print "xn   yn     e_abs  e_rel"
for x in arange(x0,xn+(h/2),h):
        e_rel = abs(y - fn_exact(x))
	print "%.2lf" % x, "%.4lf" % y, "%.4lf" % e_rel, "%lf" % (100 * (e_rel / fn_exact(x)))
	y = y + h * fn(x, y)

xn yn e_abs e_rel 1.00 1.0000 0.0000 0.000000 1.10 1.0200 0.0012 0.119666 1.20 1.0424 0.0025 0.243292 1.30 1.0675 0.0040 0.371245 1.40 1.0952 0.0055 0.503887 1.50 1.1259 0.0073 0.641577 1.60 1.1597 0.0092 0.784669 1.70 1.1968 0.0113 0.933513 1.80 1.2375 0.0136 1.088453 1.90 1.2820 0.0162 1.249829 2.00 1.3307 0.0191 1.417977 2.10 1.3839 0.0224 1.593227 2.20 1.4421 0.0261 1.775903 2.30 1.5055 0.0302 1.966323 2.40 1.5748 0.0348 2.164800 2.50 1.6504 0.0401 2.371641 2.60 1.7329 0.0460 2.587145 2.70 1.8230 0.0527 2.811607 2.80 1.9214 0.0604 3.045311 2.90 2.0290 0.0690 3.288538 3.00 2.1467 0.0788 3.541558 3.10 2.2755 0.0900 3.804635 3.20 2.4166 0.1027 4.078026 3.30 2.5713 0.1173 4.361976 3.40 2.7410 0.1339 4.656724 3.50 2.9274 0.1529 4.962500 3.60 3.1323 0.1746 5.279524 3.70 3.3578 0.1995 5.608006 3.80 3.6063 0.2281 5.948147 3.90 3.8804 0.2609 6.300138 4.00 4.1830 0.2987 6.664158 4.10 4.5177 0.3421 7.040377 4.20 4.8881 0.3923 7.428954 4.30 5.2987 0.4501 7.830035 4.40 5.7544 0.5170 8.243756 4.50 6.2608 0.5944 8.670240 4.60 6.8243 0.6840 9.109600 4.70 7.4521 0.7879 9.561934 4.80 8.1526 0.9086 10.027329 4.90 8.9352 1.0489 10.505859 5.00 9.8109 1.2123 10.997586 5.10 10.7920 1.4027 11.502555 5.20 11.8928 1.6249 12.020803 5.30 13.1296 1.8846 12.552349 5.40 14.5214 2.1885 13.097201 5.50 16.0897 2.5446 13.655352