function out = q(x,y)
      d1 = sqrt( (x+1).^2 + y.^2 );
      d2 = sqrt( (x-1).^2 + y.^2);
      l = .5; k1 = 1; k2 = 2; 
      m = 1;

      out = k1*(d1-l).*y./d1 +k2*(d2-l).*y./d2 + 9.8*m;