REM gb1bbc.bas - jackord@kw.igs.net - revised 27 Feb 08 - BBC BASIC (RTR)
      REM fieldlines and equipotentials for an equilateral -Q, +2Q, -Q array
      
      REM Use 64-bit floats
      *FLOAT 64
      
      REM Install libraries
      INSTALL @lib$+"WINLIB5"
      
      REM Initialize Window
      WindowWidth=408    : REM pixel scale 0-408
      WindowHeight=336   : REM pixel scale 0-336
      VDU 23,22,WindowWidth;WindowHeight;8,15,16,128
      button% = FN_button("Plot", 6, 6, 40, 18, FN_setproc(PROCplot), 0)
      SYS "SetWindowText", @hwnd%, "Fieldlines and Equipotentials"
      
      REM Initialize and Plot Equilateral Array
      DIM p(4): DIM q(4): DIM m(4): DIM r(4)
      d=1
      p(1)=121: q(1)=120: p(2)=204: q(2)=264: p(3)=287: q(3)=120
      m(1)=0-1: m(2)=2: m(3)=m(1)
      FOR i=1 TO 3 STEP 2
        GCOL 0 : CIRCLE p(i)*2, q(i)*2, 24
        GCOL 13 : FILL p(i)*2, q(i)*2
        GCOL 0 : VDU 5 : MOVE p(i)*2-16, q(i)*2+16 : PRINT "-Q";
      NEXT i
      GCOL 0 : CIRCLE p(2)*2, q(2)*2, 24
      GCOL 7 : FILL p(2)*2, q(2)*2
      GCOL 0 : VDU 5 : MOVE p(2)*2-14, q(2)*2+16 : PRINT "2Q";
      
      REPEAT
        WAIT 1
      UNTIL FALSE
      
      DEF PROCplot
      REM Plot Fieldlines
      GCOL 12
      FOR i=1 TO 24
        t=(2*i-1)*PI/24: x=p(2)+12*COS(t): y=q(2)+12*SIN(t)
        PROCfdln
      NEXT i
      
      REM Plot Equipotentials
      GCOL 9
      FOR i=0-11 TO 11 STEP 2
        x=204: y=168+6*i
        PROCeqpt
      NEXT i
      ENDPROC
      
      DEF PROCfdln : REM Fieldline Subroutine
      MOVE x*2, y*2
      dx=0: dy=0: PROCgrav
      dx=0-gx*d: dy=0-gy*d
      WHILE r(1)>1728 AND r(3)>1728
        PROCgrav
        dx=0-gx*d: dy=0-gy*d: x=x+dx: y=y+dy
        DRAW x*2, y*2
      ENDWHILE
      ENDPROC
      
      DEF PROCeqpt : REM Equipotential Subroutine
      MOVE x*2, y*2
      dx=0: dy=0: PROCgrav
      dx=gy*d: dy=0-gx*d: dxi=dx: xold=x
      WHILE (xold-204)*dxi>=0 OR (x-204)*dxi<0
        PROCgrav
        dx=gy*d: dy=0-gx*d: xold=x: x=x+dx: y=y+dy
        DRAW x*2, y*2
      ENDWHILE
      ENDPROC
      
      DEF PROCgrav : REM Field Subroutine
      xx=x+dx/2: yy=y+dy/2: gx=0: gy=0
      FOR k=1 TO 3
        r(k)=((xx-p(k))*(xx-p(k))+(yy-q(k))*(yy-q(k)))^1.5
        gx=gx+m(k)*(p(k)-xx)/r(k): gy=gy+m(k)*(q(k)-yy)/r(k)
      NEXT k
      gg=(gx*gx+gy*gy)^.5: gx=gx/gg: gy=gy/gg
      ENDPROC