PROGRAM rc
! simulation of RC circuit with ac voltage source
LIBRARY "csgraphics"
CALL initial(Q,R,tau,V0,omega,tmax,t,dt)
CALL set_up_windows(#1,#2,V0,tmax)
DO while t <= tmax
   CALL scope(Q,I,R,tau,V0,omega,t,dt)
   LET t = t + dt
   CALL source_voltage(#1,V0,omega,t)
   CALL output_voltage(#2,I,R,t)
LOOP
CLOSE #1
CLOSE #2
END

SUB initial(Q,R,tau,V0,omega,tmax,t,dt)
    LET t = 0
    LET V0 = 1                    ! amplitude of external voltage
    INPUT prompt "external voltage frequency (Hertz) = ":f
    LET omega = 2*pi*f            ! angular frequency
    INPUT prompt "resistance (ohms) = ": R
    INPUT prompt "capacitance (farads) = ": C
    INPUT prompt "time step dt = ": dt
    LET Q = 0
    LET tau = R*C                 ! relaxation time
    LET period = 1/f              ! period of external frequency
    IF period > tau then
       LET tmax = 2*period
    ELSE
       LET tmax = 2*tau
    END IF
END SUB

SUB set_up_windows(#1,#2,V0,tmax)
    OPEN #1: screen 0,1,0,0.5
    LET tmin = 0
    LET Vmin = -V0
    CALL draw_axes(tmin,tmax,Vmin,V0)
    PRINT "source voltage"
    OPEN #2: screen 0,1,.5,1
    CALL draw_axes(tmin,tmax,Vmin,V0)
    PRINT "voltage drop across resistor"
END SUB

SUB scope(Q,I,R,tau,V0,omega,t,dt)
    ! compute voltage drops
    DECLARE DEF V
    LET I = V(V0,omega,t)/R - Q/tau
    LET Q = Q + I*dt
END SUB

DEF V(V0,omega,t) = V0*cos(omega*t)

SUB source_voltage(#1,V0,omega,t)
    DECLARE DEF V
    WINDOW #1
    SET COLOR "blue"
    PLOT LINES: t,V(V0,omega,t);
END SUB

SUB output_voltage(#2,I,R,t)      ! voltage drop across the resistor
    WINDOW #2
    SET COLOR "red"
    PLOT LINES: t,I*R;
END SUB