Created
February 14, 2019 15:53
-
-
Save csvance/02cb7f50008f26947a236c93eb1b26a0 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
clf | |
hold on | |
disp("Note: You must define the CX output as the general solution here after running one time for this code to work") | |
syms y(x) x d C2 C4 | |
xrange = -5:0.1:5; | |
yrange = -10:0.1:10; | |
% Quiver | |
%ode = @(y, x) (x+2*y^2)/(2*y*x); | |
ode = @(y, x) (y+x)/(y-x); | |
odev = vectorize(ode); | |
odev = eval(odev); | |
[X,Y] = meshgrid(xrange, yrange); | |
S = feval(odev, Y, X); | |
L = sqrt(1+S.^2); | |
quiver(X, Y, 1./L, S./L, 0.5); | |
% General Solution | |
%solutions = dsolve(diff(y)==(x+2*y^2)/(2*y*x)); | |
solutions = dsolve(diff(y)==(y+x)/(y-x)); | |
% After dsolve, we want y to be not be y(x) | |
syms y | |
disp(solutions); | |
% IVPs | |
ivp_problem = true; | |
ivp_y = [-2]; | |
ivp_x = [-2]; | |
% If not normal IVP use values of C | |
Cs = [0 1 2 3]; | |
if ivp_problem | |
sol_count = length(ivp_y); | |
else | |
sol_count = length(Cs); | |
end | |
leg = {'Direction'}; | |
plot_idx = 2; | |
last_sol = false; | |
for idx = 1:length(solutions) | |
eq = y == solutions(idx); | |
sol = solve(eq, C4); | |
if sol == last_sol | |
continue; | |
end | |
last_sol = sol; | |
disp(sol) | |
for ivp = 1:sol_count | |
if ivp_problem | |
C = eval(subs(sol, {y x}, {ivp_y(ivp), ivp_x(ivp)})); | |
else | |
C = Cs(ivp); | |
end | |
color = rand(1,3); | |
disp(sol); | |
fcontour(sol, [xrange(1) xrange(length(xrange)) yrange(1) yrange(length(yrange))],'LineWidth',2,'LineColor',color,'LevelList', [C]); | |
if ivp_problem | |
leg{plot_idx} = sprintf("%s: y(%.2f) == %.2f", char(sol), ivp_y(ivp), ivp_x(ivp)); | |
else | |
leg{plot_idx} = sprintf("%s: C == %.2f",char(sol), C); | |
end | |
plot_idx = plot_idx + 1; | |
end | |
end | |
legend(leg, 'FontSize', 10) | |
title(sprintf("ODE: %s", char(ode))); | |
axis([xrange(1) xrange(length(xrange)) yrange(1) yrange(length(yrange))]) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment