diff --git a/rk3.py b/rk3.py new file mode 100644 index 0000000..5b6201a --- /dev/null +++ b/rk3.py @@ -0,0 +1,30 @@ +#!/usr/bin/env python +import math + +def rk3(h,tn,yn,ys,param=[]): + k1 = ys(tn,yn,*param) + k2 = ys(tn+h/2,yn+h*k1/2,*param) + k3 = ys(tn+h,yn - h*k1 + 2*h*k2, *param) + return tn+h, yn + h*(k1/3 + k2/3 + k3/3) + +def main(): + U = lambda tn,f : math.sin(2*math.pi*tn*f) + R = 5e3 + C = 1e-6 + us = lambda tn, un, f: (U(tn,f) - un)/(R*C) + fc = 1.0/(2*math.pi*R*C) + for f in [ 2**(i) * fc for i in range(-5,8)]: + ucmax = 0 + ucmin = 0 + tn = 0 + un = 0 + for i in range(1000000): + tn, un = rk3(0.000001, tn, un, us, [f]) + ucmax = max(un,ucmax) + ucmin = min(un,ucmin) + + print(" %10.3f : %10.6f %10.6f" % (f,ucmax,ucmin)) + + +if __name__ == '__main__': + main()