Plotting a 3d figure for a function with 4 variables?

Hi, I have a function with 4 variables (x,y,z,w) which means d(x,y,z,w)=0, and I have to plot a 3D figure for just 3 variable (x,y,z) without specific values. could you please tell me which matlab function is useful?

here is my code:

function [] = critical()
[x,y,z] = meshgrid(-4:0.1:4,-4:0.1:4,-4:0.1:4);
w=1;
%parameter
a=25;K=0.5;K_p=0.4;K_b=0.3; V_s=0.8;K_s=0.2; gamma=6;K_bar=1.5e-4;

%Function
phi_1 = y.^3 ./ ( y.^3 + K ^3 );
phi_2 = w.^2 ./ ( K_p^2 + w.^2);
phi_3= K_p^2 ./ ( K_p^2 + w.^2);
teta = phi_1 .* phi_2 .* z;
aa = phi_3 .* ( 1 - phi_3).*x;
b = teta./(teta+K_b*(teta+aa));
Beta=(1-K_bar*(gamma*(x-y)./y).^2)./(1+(K_s./y).^2);
alpha=a*b*gamma.*(x-(1+(1/gamma).*y))/V_s;

d=alpha-Beta;
fv=isosurface(x,y,z,d,0)
fv = isosurface(d,0)
fvc = isosurface(d,'blue')
pp = patch(isosurface(x,y,z,d,-0))
isonormals(x,y,z,d,pp)
set(pp,'FaceColor','red','EdgeColor','none');
view(3); 
axis tight
camlight

end

when I define w as a vector [w]=meshgrid(-4:0.1:4); it has this error:

Error using .*

Matrix dimensions must agree.

Error in criticalaki (line 11)

teta = phi_1 .* phi_2 .* z;

so I put w=1 which is not true and give me the figure below. but I want a soft figure like surface.

[darova]

Are your x,y,z matrices 3D?

[Adam Danz]

> ... but it doesn't give me a nice figure

What isn't nice about it? How does it not address your goal? Is the problem associated with the creation of patch values (i.e. your use of isosurface) or is the problem with plotting the values? We don't have sufficient information to address the question and the question itself is fuzzy.

Clearly state the goal, the problem you're having achieving the goal, and supply any data and code needed to reproduce the problem if necessary.

[M]

Thanks for your comments, I edited my question.

[darova]

Do you have original formulas for w?

[M]

No, it is like other variables, doesn't have any formula ...

[Adam Danz]

Your code is copied below so that Matlab Central's Run Feature can produce the plot. I also added 4 semicolons to suppress unnecessary outputs. Lastly, I added the line you described to define w with meshgrid() but I do not get that error.

[x,y,z] = meshgrid(-4:0.1:4,-4:0.1:4,-4:0.1:4);
% w=1;
[w]=meshgrid(-4:0.1:4)

w = 81×81

-4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000 -4.0000 -3.9000 -3.8000 -3.7000 -3.6000 -3.5000 -3.4000 -3.3000 -3.2000 -3.1000 -3.0000 -2.9000 -2.8000 -2.7000 -2.6000 -2.5000 -2.4000 -2.3000 -2.2000 -2.1000 -2.0000 -1.9000 -1.8000 -1.7000 -1.6000 -1.5000 -1.4000 -1.3000 -1.2000 -1.1000

%parameter
a=25;K=0.5;K_p=0.4;K_b=0.3; V_s=0.8;K_s=0.2; gamma=6;K_bar=1.5e-4;
%Function
phi_1 = y.^3 ./ ( y.^3 + K ^3 );
phi_2 = w.^2 ./ ( K_p^2 + w.^2);
phi_3= K_p^2 ./ ( K_p^2 + w.^2);
teta = phi_1 .* phi_2 .* z;
aa = phi_3 .* ( 1 - phi_3).*x;
b = teta./(teta+K_b*(teta+aa));
Beta=(1-K_bar*(gamma*(x-y)./y).^2)./(1+(K_s./y).^2);
alpha=a*b*gamma.*(x-(1+(1/gamma).*y))/V_s;
d=alpha-Beta;
fv=isosurface(x,y,z,d,0);
fv = isosurface(d,0);
fvc = isosurface(d,'blue');
pp = patch(isosurface(x,y,z,d,-0));
isonormals(x,y,z,d,pp)
set(pp,'FaceColor','red','EdgeColor','none');
view(3); 
axis tight
camlight

[M]

I copied your code but it still gives me an error...

Anyway thanks for your time, Adam.

[Adam Danz]

What Matlab release are you using?

phi_1 and z should be size 81x81x81, phi_2 sould be 81x81.

[M]

R2015b(8.6.0.267246)

Adam, do you know other matlab function to use instead of "isosurface" for this problem?

I mean first we use S = solve(d,y) and then plot3d figure for this, so what is your suggestion?

[Adam Danz]

I see. The problem must be due to your older release, then. Explicit expansion of variables is supported from R2016b onward.

I am not familiar enough with your data to recommend an alternative to isosurface.

[DGM]

This should at least fix your errors.

%teta = phi_1 .* phi_2 .* z;
teta = phi_1 .* bsxfun(@times,phi_2,z);
%aa = phi_3 .* ( 1 - phi_3).*x;
aa = bsxfun(@times,phi_3 .* ( 1 - phi_3),x);

isosurface doesn't tend to give a very pretty figure for a coarse grid. In order to use surf(), you'd need to express (either literally or mathematically) the location of the isosurface as a mesh of points in x,y,z space. It might look a lot nicer than using isosurface(), but it's often a lot more complicated to accomplish.

FWIW, you can always throw extra memory at it and use a finer mesh. I don't have much memory, so I just did this:

vv = linspace(-4,4,400);
[x,y,z] = meshgrid(single(vv));
[w]=meshgrid(vv);

%parameter
a=25;K=0.5;K_p=0.4;K_b=0.3; V_s=0.8;K_s=0.2; gamma=6;K_bar=1.5e-4;

%Function
phi_1 = y.^3 ./ ( y.^3 + K ^3 );
phi_2 = w.^2 ./ ( K_p^2 + w.^2);
phi_3= K_p^2 ./ ( K_p^2 + w.^2);
teta = phi_1 .* bsxfun(@times,phi_2,z);
aa = bsxfun(@times,phi_3 .* ( 1 - phi_3),x);
clear phi_1 phi_2 phi_3 w
b = teta./(teta+K_b*(teta+aa));
Beta=(1-K_bar*(gamma*(x-y)./y).^2)./(1+(K_s./y).^2);
alpha=a*b*gamma.*(x-(1+(1/gamma).*y))/V_s;
d=alpha-Beta;
clear alpha Beta a b aa teta

pp = patch(isosurface(x,y,z,d,-0));
isonormals(x,y,z,d,pp)
set(pp,'FaceColor','red','EdgeColor','none');
view(3); 
axis tight
camlight
grid on

It's a bit better, but it's still rough. At this resolution, each one of the 3D arrays occupies about 250MB, so it gets impractical to do much more.

[Adam Danz]

+1 DGM

Move to answer?

[DGM]

I kind of assumed that the core of the question was "how to I rewrite this to make a surf plot", and I really left that unanswered.

[M]

Thanks Adam and DGM for your comments,

DGM, I copied your code and run it, it has an error that "out of memory"...

I did something else, solve "alpha-Beta=0" and found an equation that x is written based on y,z,w and use "surface" or "plot3" to plot it.

[DGM]

If you're running out of memory, with that example, then you'd need to reduce the number of points to something less than 400.

If you managed to get it working with surf, then you probably don't need to bother with my example.