I am trying to graph the direction field of dy/dx=-x/y and the solution curve that satisfies the initial conditions y(0)=4. The solution curve doesn't go left of the y-axis.

I tried changing the 0 in the ts=[0 maxx] to a negative value and that results in the code not running forever without graphing anything.

1;

% Preparatory Code
clear all;
hold on
grid minor

% Set the limits of the graph

minx = -3.;
maxx = 5.;
miny = -3.;
maxy = 5.;

% Step size
h = 0.5;

% Make a meshgrid for the direction field

[x, y]=meshgrid(minx:h:maxx, miny:h:maxy);

% Define the differential equation
dt = ones(size(x));
dp = -(x./y);

% Find the angle of the diff eq
angle = atan(dp);
xL = cos(angle);
yL = sin(angle);

h2 = quiver(x, y, xL, yL);


set(h2, 'AutoScale', 'on', 'AutoScaleFactor', 0.5, 'ShowArrowHead', 'off');
axis ([minx, maxx, miny, maxy], 'square');

% Initial Value 4
ts = [0 maxx];
y0 = [4; 0];
ode = @(ts,p) firstO(ts,p);
[ts,y] = ode45(ode, ts, y0);
plot(ts,y(:,1),"Color",[1 0 0]);

% Axis Labels
xlabel('X')
ylabel('Y')
title('Direction field of dy/dx=-x/y')


hold off;

function dpdt = firstO(x,y)
dpdt = -(x./y);
end

[James Tursa]

I'm assuming the dx/dy = -x/y is a typo and you really meant dy/dx = -x/y?

[Jared Shepard]

Yes, I meant dy/dx = -x/y.

Answer 1

[Gyan Vaibhav]

Hi Jared,

In the above code, there were issues with the use of the “ode45” function. This function is used for solving initial value problems of the form dy/dt = f(t,y) with a given initial value y(t0) = y0. However, the differential equation in the question is in the form dy/dx = -x/y, which is not in the standard form for “ode45”.

In addition, the initial condition provided was y(0) = 4, but in the code, a 2-element vector [4; 0] was given as the initial condition, which was incorrect.

The differential equation dy/dx = -x/y with the initial condition y(0) = 4 can be solved using the “dsolve” function in MATLAB. However, due to the singularity at y = 0, “dsolve” returns a complex solution. When plotting the real part of this solution, it ends up being a straight line at y = 0.

Here is the MATLAB code for it:

% Preparatory Code
clear all;
hold on
grid minor
% Set the limits of the graph
minx = -3.;
maxx = 5.;
miny = -3.;
maxy = 5.;
% Step size
h = 0.5;
% Make a meshgrid for the direction field
[x, y] = meshgrid(minx:h:maxx, miny:h:maxy);
% Define the differential equation
dp = -(x./y);
% Find the angle of the diff eq
angle = atan(dp);
xL = cos(angle);
yL = sin(angle);
h2 = quiver(x, y, xL, yL);
set(h2, 'AutoScale', 'on', 'AutoScaleFactor', 0.5, 'ShowArrowHead', 'off');
axis ([minx, maxx, miny, maxy], 'square');
% Solve differential equation using dsolve
syms x y(x)
Dy = diff(y, x);
eqn = Dy == -x/y;
cond = y(0) == 4;
soln = dsolve(eqn, cond);
% Plot the real part of the solution
fplot(real(soln), [0, maxx], 'Color', [1 0 0]);
% Axis Labels
xlabel('X')
ylabel('Y')
title('Direction field of dy/dx=-x/y')
hold off;

Alternatively, you can directly solve this differential equation by recognizing that it's a first-order homogeneous equation. Its general solution is y^2 + x^2 = C. Given the initial condition y(0) = 4, we find C = 16. So, the particular solution that satisfies the initial condition is y = sqrt(16 - x^2), for x >= 0.

Here is the code with this approach:

clear all;
hold on
grid minor
% Set the limits of the graph
minx = -3.;
maxx = 5.;
miny = -3.;
maxy = 5.;
% Step size
h = 0.5;
% Make a meshgrid for the direction field
[x, y] = meshgrid(minx:h:maxx, miny:h:maxy);
% Define the differential equation
dp = -(x./y);
% Find the angle of the diff eq
angle = atan(dp);
xL = cos(angle);
yL = sin(angle);
h2 = quiver(x, y, xL, yL);
set(h2, 'AutoScale', 'on', 'AutoScaleFactor', 0.5, 'ShowArrowHead', 'off');
axis ([minx, maxx, miny, maxy], 'square');
% Solution curve that satisfies the initial condition y(0)=4
x_sol = 0:h:maxx;
y_sol = sqrt(16 - x_sol.^2);
plot(x_sol, y_sol, "Color", [1 0 0]);
% Axis Labels
xlabel('X')
ylabel('Y')
title('Direction field of dy/dx=-x/y')
hold off;

Though, the differential equation dy/dx = -x/y is not in the standard form for “ode45”. The “ode45” function in MATLAB is used for solving initial value problems of the form dy/dt = f(t,y) with a given initial value y(t0) = y0.

However, you can rewrite your equation in a form that “ode45” can handle by introducing a new variable. Let's say z = y/x, then y = zx. Differentiating both sides with respect to x gives dy/dx = z + x*dz/dx.

Now you can substitute this into your original equation to get a new equation in terms of z and x:

z + x*dz/dx = -x/(zx) dz/dx = -1/x^2 - z/x

You can refer to the documentation links regarding more information on “ode45” and “dsolve” function:

• https://mathworks.com/help/matlab/ref/ode45.html
• https://mathworks.com/help/symbolic/dsolve.html

I hope these suggestions help you resolve the issue. 

Best regards, 

Gyan