function LienardWiechertFastElectron
eps0=8.854e-12; qe=-1.6e-19; c=299.792e6;
v0=0.91*c; % Velocity of the electron

% Meshgrid
T=1e-4; 
dx=c*T/50;  % Distance of the points x-direction
dy=c*T/50;
x=-c*T:dx:c*T;
y=-c*T:dy:c*T;
[X,Y]=meshgrid(x,y); %

    for i=1:length(x)
        for j=1:length(y)
            r=[X(i,j) Y(i,j)];
            phi(i,j)=LWP(r,qe,-1/4/pi/eps0,c,T/2);
        end
    end
    
    % Drawing the figure
    colormap jet
    phimax=-1e-12;
    surfc(X/1000,Y/1000,max(phi,phimax));
    grid on
    xlabel('x [km]','FontWeight','bold','FontSize',13)
    ylabel('y [km]','FontWeight','bold','FontSize',13)
    zlabel('\phi [J/C]','FontWeight','bold','FontSize',13)
    view([-11.7 15.5789781021898]);


    function phi=LWP(r,ms,gamma,c,t)
        rsr=location(t); 
        dt=sqrt((rsr-r)*(rsr-r)')/c;
        for ii=1:1000 % Iteration for retarded time
            dtn=dt;
            rsr=location(t-dt);
            dr=r-rsr;
            absdr=sqrt(dr*dr');
            dt=absdr/c;
            if abs(dt-dtn)<1e-8
                break
            end
        end
        beta=velocity(t-dt)/c;
        phi=-gamma*ms/(absdr-dr*beta');
    end

    function rs=location(t)
        rs=[v0 0]*t;
    end
    function vs=velocity(t)
        vs=[v0 0];
    end
end