function ft = f(at)
    global kvec filt epsilon
    u=real(fft(at.*filt));
    u2=u.*u/2;
    b=ifft(u2);
    ft=-kvec.*kvec*epsilon.*at+1i*kvec.*b;
end
