MATLAB Answers

How to plot a multivariable function with constraints

4 views (last 30 days)
Lorenzo Manni
Lorenzo Manni on 28 Aug 2021
Commented: Lorenzo Manni on 30 Aug 2021
Hello,
I plotted a 2-variables function and I would like to highlight (or plot/intersect) only the z points for which x (k in the function) and y satisfy the following constraints: 250020*((pi/4)*k^2*y)-325.33=0
I have been trying several approaches but I cannot make it happen. I hope you can help.
Thank you!
f=@(k,y) (8120*(((pi*((((k+(2*(5*0.0026))+(2*0.0026))+(((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350))*2)+((4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*2)).^2)-((((k+(2*(5*0.0026))+(2*0.0026))+(((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350))*2)+((4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*2))-(2*((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350)))).^2))*y)/4)+(((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))*(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*y.*24)))+(8120*(pi*((k.^2)-((k-(2*((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350)))).^2))*y)/4)+(7450*((pi*y.*(((k+(2*(5*0.0026))).^2)-(k.^2)))/4)+(8933*pi*(0.003/2).^2*((15*(2*4))*((2*y)+((2*(pi*((k+(2*(5*0.0026))+(2*0.0026))+(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350)))))))/8)+(4*(((pi*((k+(2*(5*0.0026))+(2*0.0026))+(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))))*3));
fsurf(f);
xlim([0.02 0.1]);
ylim([0.04 0.4]);
zlim([0 100]);

Answers (1)

John D'Errico
John D'Errico on 28 Aug 2021
Can you solve for y, as a function of k, given the constraint? (YES.)
If so, then do you really have a TWO variable problem? (NO.)
syms k y
ysol = solve(250020*((pi/4)*k^2*y)-325.33==0,y)
ysol = 
z = (8120*(((pi*((((k+(2*(5*0.0026))+(2*0.0026))+(((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350))*2)+((4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*2)).^2)-((((k+(2*(5*0.0026))+(2*0.0026))+(((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350))*2)+((4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*2))-(2*((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350)))).^2))*y)/4)+(((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))*(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*y.*24)))+(8120*(pi*((k.^2)-((k-(2*((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350)))).^2))*y)/4)+(7450*((pi*y.*(((k+(2*(5*0.0026))).^2)-(k.^2)))/4)+(8933*pi*(0.003/2).^2*((15*(2*4))*((2*y)+((2*(pi*((k+(2*(5*0.0026))+(2*0.0026))+(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350)))))))/8)+(4*(((pi*((k+(2*(5*0.0026))+(2*0.0026))+(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))))*3));
z_k = subs(z,y,ysol)
z_k = 
It may look messy, but who cares?
fplot(z_k,[0.02,0.1])
xlabel 'k'
ylabel 'z'
Again, this is not a 2 variable problem, due to the constraint.
Can you plot it in 3-dimensions anyway? Well, yes.
K = linspace(0.02,0.1,100);
yfun = matlabFunction(ysol);
Y = yfun(K);
f = @(k,y) (8120*(((pi*((((k+(2*(5*0.0026))+(2*0.0026))+(((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350))*2)+((4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*2)).^2)-((((k+(2*(5*0.0026))+(2*0.0026))+(((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350))*2)+((4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*2))-(2*((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350)))).^2))*y)/4)+(((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))*(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))*y.*24)))+(8120*(pi*((k.^2)-((k-(2*((pi*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(8*0.97*2.350)))).^2))*y)/4)+(7450*((pi*y.*(((k+(2*(5*0.0026))).^2)-(k.^2)))/4)+(8933*pi*(0.003/2).^2*((15*(2*4))*((2*y)+((2*(pi*((k+(2*(5*0.0026))+(2*0.0026))+(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350)))))))/8)+(4*(((pi*((k+(2*(5*0.0026))+(2*0.0026))+(4.5*(((pi*(k+(2*(5*0.0026))+(2*0.0026)))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))))/24)-((2*pi()*((k+(2*(5*0.0026)))/2)*((1.49*(5*0.0026))/((5*0.0026)+0.0026)))/(24*0.97*2.350))))))*3));
f = str2func(vectorize(func2str(f)));
plot3(K,Y,f(K,Y),'o-')
grid on
box on
xlabel 'k'
ylabel 'y'
zlabel 'z'
  4 Comments
Lorenzo Manni
Lorenzo Manni on 30 Aug 2021
The issue is that now the line
f = str2func(vectorize(func2str(f)))
does not work because there is no
f=@(k,y)
before, since I am writing the function using
syms k y
How would you fix that @John D'Errico?
Many thanks

Sign in to comment.

Products


Release

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!