Is there a way to vpasolve elementwise?

Hello,

I am hoping someone can point me in the right direction with regards to the following problem. I have an equation with two variables, x and y, and an unknown c, and I'd like to solve the equation for c over a range of x and y values. The equation cannot be solved for c analytically. Eventually, I'd like to plot the solutions for c as function of x and y in a contour plot. Below is the code that I thought would work for a much simplified example of the actual problem; this may help define the problem.

x = 1:5;
y = 1:5;
[X,Y] = meshgrid(x,y);
syms c
sol = vpasolve(X + c*Y - 3 == 0,c);

Obviously, this did not work. I'm sure I'm thinking about this wrong way. Any helpful ideas are welcome.

Thank you, John.

Answer 1

[Star Strider]

The contour plot is actually the solution. Solving a two-variable problem like that is one of its hidden benefits.

I’m not quite sure what you want to do, but I solved for ‘c’ here, and then had contour provide the (x,y) coordinates for that over the grid you gave it. Here, I asked it for the coordinates corresponding to c=0:

x = 1:5;
y = 1:5;
[X,Y] = meshgrid(x,y);
F = (3-X)./Y;
sol = contour(X, Y, F, [0 0]);
xy = sol(:,2:end)

The ‘xy’ matrix are the coordinates, x = xy(1,:) and y = xy(2,:). I refer you to the documentation on contour for the details. Substitute other values for the fourth argument to define other values of ‘c’.

[John Moseley]

Thank you very much for your fast response. This definitely helps. My issue is that in the actual equation I'm using with x, y, and c, there is not an analytical solution for c. So, I was trying to use vpasolve.

John

[Star Strider]

My pleasure.

I figured that this must be a proxy for a different problem. If you have one equation in two unknowns, while you cannot have a unique solution for ‘c’, you can obviously get a range of solutions for ‘c’ if you define ranges for the other two variables.

A problem arises if powers of ‘c’ are involved, because there will be as many values of ‘c’ for every (X,Y) pair as there are powers of ‘c’, and they may not all be real.

[John Moseley]

For a given x,y pair, there are three solutions for c. There is only one positive solution, and that's the one I need, so I can write that in the code easy enough. My trouble is generating all the positive solutions for c over the [X,Y] mesh.

[Star Strider]

Without knowing your actual equation, the only other option I can offer is to loop over the individual elements in the mesh:

fn = @(c,X,Y) X + c*Y - 3;
x = 1:5;
y = 1:5;

for k1 = 1:length(x)
    for k2 = 1:length(y)
        C(k1,k2) = fsolve(@(c) fn(c,x(k1),y(k2)), 1);
    end
end

You can also use fzero instead of fsolve.

[John Moseley]

Perfect, Thanks!

[Star Strider]

My pleasure!