How can I create a 2D contour Plot from X,Y,Z table?

I am trying to plot the 2D contour from the equation. I have tried my code. However, It does not work. I give my code below. Thank you so much!.

     num=50;
    Lpmax=700e-6; 
    Lpmin=20e-6; 
    Lr=2.77e-3;
    Cr=2.5e-9;
    
    f=1/(2*pi*sqrt(Lr*Cr*0.5));
    w=2*pi*f;
    T=1/f;
    Tdead=T*0.011;    
    Lp_x11=linspace(Lpmin,Lpmax,num);
    Vs=linspace(140,340,num);
    Vo=180;
    Re_x11=200;
    theta=0;
    for i=1:num
    % solve the D
    syms D11
    X=solve(sin(D11*pi)/(2*(1-D11))-Vo/Vs(i));
    D_ch=X;
    %
    theta_dch=(2*pi*Tdead)/T;
    A2=(1./(Re_x11.*(pi.^2).*(1-D_ch))).*2.*w.*Lp_x11(i).*sin(pi.*D_ch).*sin(pi.*(1-D_ch));
    Qlr1=(2.*Vs(i).*sin(pi.*D_ch)./(Re_x11.*pi.*(1-D_ch).*w)).*(sin(-pi.*D_ch-theta_dch)-sin(-pi.*D_ch));
    Qlr2=(2.*Vs(i).*sin(pi.*D_ch)./(Re_x11.*pi.*(1-D_ch).*w)).*(sin(pi*D_ch-theta_dch)-sin(pi.*D_ch));
    Qlp_S1=-(1./w).*(Vs(i).*theta_dch./(w.*Lp_x11(i))).*(A2./(1-D_ch)-pi.*D_ch+theta_dch./2);
    Qlp_S2=(1./w).*(D_ch.*Vs(i).*theta_dch)./(w.*Lp_x11(i).*(1-D_ch)).*(A2./D_ch+pi.*(1-D_ch)-theta_dch./2);
    Qs1(i)=Qlp_S1+Qlr1;
    Qs2(i)=Qlr2+Qlp_S2;
    end 
    xv = linspace(min(Lp_x11), max(Lp_x11), 1000);
    zv = linspace (min (Vs), max (Vs), 1000);
    [X,Z] = meshgrid(xv, zv);
    P = griddata(Lp_x11, Vs, Qs1, X, Z);
    
    figure
    scatter3(x, z, p, '.')
    grid on
    figure
    meshc (X, Z, P)
    grid on
    figure
    contourf(X, Z, P)
    grid on

Answer 1

[Cris LaPierre]

The first issue is that one of your inputs to griddata is symbolic. It must be numeric. Use the double function to convert Qs1 to doubles.

P = griddata(Lp_x11, Vs, double(Qs1), X, Z);

Next, your scatter3 function call does not use variables that exist. Variables are case sensitive. Use X, Z and P. In addition, the inputs to scatter3 must be vectors, not arrays. Use a colon to turn your arrays into vectors.

scatter3(X(:), Z(:), P(:), '.')

Surface plots (meshes, contours) will not work with the data as you currently have it since the data is more of a line than a surface. You have gridded the data, but most of the values are NaN. Scatter3 may be the best choice.

You could also try just plotting in 3D.

plot3(Lp_x11, Vs, double(Qs1))

Upon closer inspection, you'll see that you are ony solving you equation along the diagonal (Vs(i),Lp_x11(i)). This is why your results of all NaN. You are essentially trying to determine an entire surface from a single line. Instead, you want to user your loop to pick one value (say Vs(i)), and solve for all values of Lp_x11. If you capture the results correctly, Qs1 and Qs2 will be num x num arrays.

At that point, you want to use interp2, not gridded data, to increase the resolution of your data.

Here's your code modified.

num=50;
Lpmax=700e-6;
Lpmin=20e-6;
Lr=2.77e-3;
Cr=2.5e-9;

f=1/(2*pi*sqrt(Lr*Cr*0.5));
w=2*pi*f;
T=1/f;
Tdead=T*0.011;
Lp_x11=linspace(Lpmin,Lpmax,num);
Vs=linspace(140,340,num);
Vo=180;
Re_x11=200;
theta=0;

% solve the D
syms D11
theta_dch=(2*pi*Tdead)/T;
for r=1:length(Vs)
    X=vpasolve(sin(D11*pi)/(2*(1-D11))-Vo/Vs(r));
    D_ch=X;
    A2=(1./(Re_x11.*(pi.^2).*(1-D_ch))).*2.*w.*Lp_x11.*sin(pi.*D_ch).*sin(pi.*(1-D_ch));
    Qlr1=(2.*Vs(r).*sin(pi.*D_ch)./(Re_x11.*pi.*(1-D_ch).*w)).*(sin(-pi.*D_ch-theta_dch)-sin(-pi.*D_ch));
    Qlr2=(2.*Vs(r).*sin(pi.*D_ch)./(Re_x11.*pi.*(1-D_ch).*w)).*(sin(pi*D_ch-theta_dch)-sin(pi.*D_ch));
    Qlp_S1=-(1./w).*(Vs(r).*theta_dch./(w.*Lp_x11)).*(A2./(1-D_ch)-pi.*D_ch+theta_dch./2);
    Qlp_S2=(1./w).*(D_ch.*Vs(r).*theta_dch)./(w.*Lp_x11.*(1-D_ch)).*(A2./D_ch+pi.*(1-D_ch)-theta_dch./2);
    
    % Vs will be the y axis, which corresponds to the rows of the matrix
    % Capture the entire row of data using the indexing below.
    Qs1(r,:)=Qlp_S1+Qlr1; 
    Qs2(r,:)=Qlr2+Qlp_S2;
end
xv = linspace(min(Lp_x11), max(Lp_x11), 1000);
zv = linspace (min (Vs), max (Vs), 1000);
[X,Z] = meshgrid(xv, zv);

% P = griddata(Lp_x11, Vs, double(Qs1), X, Z);
P = interp2(Lp_x11, Vs, double(Qs1), X, Z);

figure
scatter3(X(:), Z(:), P(:), '.')
grid on

figure
meshc(X, Z, P)
grid on

figure
contourf(X, Z, P)
grid on

[Tin Truong Chanh]

Thank you so much for your help. Your method works well. my problem is solved.

Answer 2

[Walter Roberson]

    num=50;
    Lpmax=700e-6; 
    Lpmin=20e-6; 
    Lr=2.77e-3;
    Cr=2.5e-9;
    
    f=1/(2*pi*sqrt(Lr*Cr*0.5));
    w=2*pi*f;
    T=1/f;
    Tdead=T*0.011;    
    Lp_x11=linspace(Lpmin,Lpmax,num);
    Vs=linspace(140,340,num);
    Vo=180;
    Re_x11=200;
    theta=0;
    Qs1 = zeros(1,num);
    Qs2 = zeros(1,num);
    for i=1:num
    % solve the D
    syms D11
    X=vpasolve(sin(D11*pi)/(2*(1-D11))-Vo/Vs(i));
    D_ch=double(X);
    %
    theta_dch=(2*pi*Tdead)/T;
    A2=(1./(Re_x11.*(pi.^2).*(1-D_ch))).*2.*w.*Lp_x11(i).*sin(pi.*D_ch).*sin(pi.*(1-D_ch));
    Qlr1=(2.*Vs(i).*sin(pi.*D_ch)./(Re_x11.*pi.*(1-D_ch).*w)).*(sin(-pi.*D_ch-theta_dch)-sin(-pi.*D_ch));
    Qlr2=(2.*Vs(i).*sin(pi.*D_ch)./(Re_x11.*pi.*(1-D_ch).*w)).*(sin(pi*D_ch-theta_dch)-sin(pi.*D_ch));
    Qlp_S1=-(1./w).*(Vs(i).*theta_dch./(w.*Lp_x11(i))).*(A2./(1-D_ch)-pi.*D_ch+theta_dch./2);
    Qlp_S2=(1./w).*(D_ch.*Vs(i).*theta_dch)./(w.*Lp_x11(i).*(1-D_ch)).*(A2./D_ch+pi.*(1-D_ch)-theta_dch./2);
    Qs1(i)=Qlp_S1+Qlr1;
    Qs2(i)=Qlr2+Qlp_S2;
    end
    
    nnz(~isfinite(Lp_x11)), nnz(~isfinite(Vs)), nnz(~isfinite(Qs1))

ans = 0

ans = 0

ans = 0

    xv = linspace(min(Lp_x11), max(Lp_x11), 1000);
    zv = linspace(min(Vs), max(Vs), 1000);
    [X,Z] = ndgrid(xv, zv);
    F = scatteredInterpolant(Lp_x11(:), Vs(:), Qs1(:), 'linear', 'linear');
    P = F(X, Z);
    
    scatter(X(:), Z(:), 10)

    figure
    scatter3(X(:), Z(:), P(:), '.')
    grid on

    figure
    meshc(X, Z, P)
    grid on

    figure
    contourf(X, Z, P, linspace(-3e-6,4e-6,8))
    grid on

 nnz(isnan(P))

ans = 999079

 surf(X, Z, P, 'edgecolor', 'none')

[Walter Roberson]

Okay, so what is happening is that when you calculate Qs1, you index Lp_x11(i) and Vs(i) . Both of those were created by linspace() with the same number of points, so they are in relationship, Vs(i) = A*Lp_x11(i) + B for some constants A and B.

Therefore, Qs1 is only being calculated along a line, and when you go to interpolate that line into 3 space, you do not have any information to interpolate off-axes, so it all interpolates as nan.

[Tin Truong Chanh]

Thank you so much for your help!