Why does my code produce a single direction vector and not a phase portrait of all the direction vectors

I was hoping to get a phase portrait but instead I get a short spiral and no other values.

% Mesh Grid in (Diamtoms(D),Zooplankton(P))-plane
[D, P] = meshgrid(-1:0.20:200, -1:.20:200);
% Parameters
K=1; s=1; e=1; r=1; H=1; d=0.75;

% System as a 2-D function
f = @(t,X) [X(1)*(r*(1-X(1)*(1/K))-(s*X(1)/1+s*H*X(1))*X(2)); 
    X(2)*(e*(s*X(1))/(1+s*H*X(1))-d)];

%Direction Filed
% Ddot is dD/dt and Pdot is dP/dt derivatives
Ddot = D-(1/K)*D*D-(s*D)/((1+s*H*D))*P;

Warning: Matrix is singular to working precision.

Pdot = P*((s*D)/(1+s*H*D))-P*d;

Warning: Matrix is singular to working precision.

% Vector Field
figure(1)
quiver(D,P,Ddot,Pdot,'k','LineWidth',1.5, 'AutoScaleFactor', 3)
hold on

Warning: Matrix is singluar to working precision

corresponding to my derivatives and my meshgrid

Answer 1

[Dyuman Joshi]

You need to use element-wise operations for defining Ddot and Pdot.

Also, starting the values from -0.8, as for D == -1, the values of Ddot and Pdot are Infinite (as 1+s*H*D == 0)

% Mesh Grid in (Diamtoms(D),Zooplankton(P))-plane
[D, P] = meshgrid(-0.8:0.20:200);
% Parameters
K=1; s=1; e=1; r=1; H=1; d=0.75;

%Direction Filed
% Ddot is dD/dt and Pdot is dP/dt derivatives
Ddot = D-(1./K).*D.*D-(s*D)./((1+s.*H.*D)).*P;
Pdot = P.*((s.*D)./(1+s.*H.*D))-P.*d;

% Vector Field
figure
quiver(D,P,Ddot,Pdot,'k','LineWidth',1.5,'AutoScaleFactor',3)