How show direction of line inside plot and ?

in my plots i want something simple to but have a good design for ploting i got some code but they not like i am looking i want something better i have a system for plotting and 2 point i want to plot it with better design and trajectory on line are shown

i want something like that but again if possible the line are shown on the line it will be so better

%% Phase Portrait Plots in MATLAB
% Creating beautiful contour plots similar to the reference images

%% Method 1: Single Plot with Your Exact Parameters
% Define the function P = (1/2)*F1*u^2 + (1/3)*F2*u^3 - (1/2)*v^2
% With F1 = 1, F2 = -1

% Define grid
u = linspace(-1.5, 2, 1000);
v = linspace(-1.5, 1.5, 1000);
[U, V] = meshgrid(u, v);

% Define parameters
F1 = 1;
F2 = -1;

% Calculate the function P
P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;

% Create figure
figure('Position', [100, 100, 800, 600]);

% Create contour plot
contour_levels = -0.4:0.1:0.4;
[C, h] = contour(U, V, P, contour_levels, 'LineWidth', 1.5);

% Enhance the plot
colormap(parula);  % You can also try: jet, hot, cool, spring, etc.
axis equal;
grid on;
grid minor;
xlabel('u', 'FontSize', 14, 'FontWeight', 'bold');
ylabel('v', 'FontSize', 14, 'FontWeight', 'bold');
title('Phase Portrait: P = \frac{1}{2}u^2 - \frac{1}{3}u^3 - \frac{1}{2}v^2', ...
      'FontSize', 16, 'FontWeight', 'bold');

% Add colorbar
colorbar('FontSize', 12);

% Set axis properties
set(gca, 'FontSize', 12, 'LineWidth', 1.2);
xlim([-1.5, 2]);
ylim([-1.5, 1.5]);

% Add equilibrium points (find where gradient = 0)
% For this system: ∂P/∂u = F1*u + F2*u^2 = u + (-1)*u^2 = u(1-u) = 0
% So u = 0 or u = 1
% ∂P/∂v = -v = 0, so v = 0
equilibrium_points_u = [0, 1];
equilibrium_points_v = [0, 0];
hold on;
plot(equilibrium_points_u, equilibrium_points_v, 'ko', 'MarkerSize', 8, ...
     'MarkerFaceColor', 'black', 'LineWidth', 2);
hold off;

%% Method 2: Multiple Subplots with Different Parameters
figure('Position', [150, 150, 1200, 400]);

% Different parameter combinations
params = [3, -3; -1, -2; 0, -2];  % [F1, F2] pairs
titles = {'(a) F_1 = 3, F_2 = -3', '(b) F_1 = -1, F_2 = -2', '(c) F_1 = 0, F_2 = -2'};

for i = 1:3
    subplot(1, 3, i);
    
    F1 = params(i, 1);
    F2 = params(i, 2);
    
    % Calculate P for current parameters
    P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;
    
    % Create contour plot
    [C, h] = contour(U, V, P, contour_levels, 'LineWidth', 1.5);
    
    % Enhance each subplot
    colormap(parula);
    axis equal;
    grid on;
    xlabel('u', 'FontSize', 12);
    ylabel('v', 'FontSize', 12);
    title(titles{i}, 'FontSize', 14, 'FontWeight', 'bold');
    
    % Set consistent axis limits
    xlim([-1.5, 2]);
    ylim([-1.5, 1.5]);
    
    % Find and plot equilibrium points for current parameters
    if F1 ~= 0
        % Solve F1*u + F2*u^2 = 0 → u(F1 + F2*u) = 0
        eq_u = [0, -F1/F2];
        eq_v = [0, 0];
        % Filter for real solutions within plot range
        valid_idx = isreal(eq_u) & eq_u >= -1.5 & eq_u <= 2;
        eq_u = eq_u(valid_idx);
        eq_v = eq_v(valid_idx);
    else
        % If F1 = 0, then F2*u^2 = 0, so u = 0
        eq_u = 0;
        eq_v = 0;
    end
    
    hold on;
    if ~isempty(eq_u)
        plot(eq_u, eq_v, 'ko', 'MarkerSize', 6, 'MarkerFaceColor', 'black');
    end
    hold off;
    
    set(gca, 'FontSize', 10);
end

%% Method 3: Enhanced Version with Better Colors and Styling
figure('Position', [200, 200, 900, 700]);

% Use a more sophisticated colormap
custom_colormap = [
    0.2, 0.1, 0.5;  % Dark blue
    0.1, 0.4, 0.8;  % Blue
    0.1, 0.8, 0.8;  % Cyan
    0.9, 0.9, 0.1;  % Yellow
    0.8, 0.4, 0.1;  % Orange
    0.5, 0.1, 0.2   % Dark red
];

% Reset to original parameters
F1 = 1;
F2 = -1;
P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;

% Create filled contour plot for better visualization
[C, h] = contourf(U, V, P, 20, 'LineStyle', 'none');
hold on;

% Add contour lines
[C2, h2] = contour(U, V, P, contour_levels, 'LineWidth', 1.2, 'LineColor', 'black');

% Apply custom colormap
colormap(custom_colormap);

% Enhanced styling
axis equal;
xlabel('u', 'FontSize', 16, 'FontWeight', 'bold');
ylabel('v', 'FontSize', 16, 'FontWeight', 'bold');
title('Enhanced Phase Portrait', 'FontSize', 18, 'FontWeight', 'bold');

% Add colorbar with custom formatting
cb = colorbar;
cb.Label.String = 'P(u,v)';
cb.Label.FontSize = 14;
cb.Label.FontWeight = 'bold';

% Plot equilibrium points
plot([0, 1], [0, 0], 'wo', 'MarkerSize', 10, 'MarkerFaceColor', 'white', ...
     'MarkerEdgeColor', 'black', 'LineWidth', 2);

% Set axis properties
set(gca, 'FontSize', 14, 'LineWidth', 1.5);
xlim([-1.5, 2]);
ylim([-1.5, 1.5]);

% Add grid
grid on;
set(gca, 'GridColor', [0.3, 0.3, 0.3], 'GridAlpha', 0.3);

hold off;

%% Method 4: Example usage of custom functions (defined at end of file)
% Uncomment the line below to use the custom function:
% create_phase_portrait(1, -1, 'Custom Phase Portrait: F_1 = 1, F_2 = -1');

%% Method 5: Phase Portrait with Direction Field and Trajectories
figure('Position', [250, 250, 1000, 800]);

% Parameters
F1 = 1;
F2 = -1;

% Define finer grid for contours
u_fine = linspace(-1.5, 2, 1000);
v_fine = linspace(-1.5, 1.5, 1000);
[U_fine, V_fine] = meshgrid(u_fine, v_fine);
P_fine = (1/2)*F1*U_fine.^2 + (1/3)*F2*U_fine.^3 - (1/2)*V_fine.^2;

% Coarser grid for vector field
u_coarse = linspace(-1.5, 2, 20);
v_coarse = linspace(-1.5, 1.5, 15);
[U_coarse, V_coarse] = meshgrid(u_coarse, v_coarse);

% Compute the vector field (gradient descent: du/dt = -∂P/∂u, dv/dt = -∂P/∂v)
% For potential systems: 
% du/dt = -∂P/∂u = -(F1*u + F2*u^2)
% dv/dt = -∂P/∂v = -(-v) = v

DU = -(F1*U_coarse + F2*U_coarse.^2);  % -∂P/∂u
DV = V_coarse;                          % -∂P/∂v

% Alternative: Hamiltonian system (symplectic)
% DU = -V_coarse;                       % ∂P/∂v  
% DV = -(F1*U_coarse + F2*U_coarse.^2); % -∂P/∂u

% Normalize vector field for better visualization
magnitude = sqrt(DU.^2 + DV.^2);
DU_norm = DU ./ (magnitude + 1e-8);
DV_norm = DV ./ (magnitude + 1e-8);

% Create filled contour plot
[C, h] = contourf(U_fine, V_fine, P_fine, 20, 'LineStyle', 'none');
hold on;

% Add contour lines
contour_levels = -0.4:0.1:0.4;
[C2, h2] = contour(U_fine, V_fine, P_fine, contour_levels, 'LineWidth', 1.2, 'LineColor', 'black');

% Add vector field
quiver(U_coarse, V_coarse, DU_norm, DV_norm, 0.5, 'Color', 'white', ...
       'LineWidth', 1.5, 'MaxHeadSize', 0.3);

% Add streamlines from various starting points
start_u = [-1.2, -0.8, -0.4, 0.4, 0.8, 1.2, 1.6];
start_v = [-1.2, -0.8, -0.4, 0, 0.4, 0.8, 1.2];
[start_U, start_V] = meshgrid(start_u, start_v);
streamline(U_coarse, V_coarse, DU, DV, start_U(:), start_V(:), [0.1, 10000]);

% Enhanced styling
colormap(parula);
axis equal;
xlabel('u', 'FontSize', 16, 'FontWeight', 'bold');
ylabel('v', 'FontSize', 16, 'FontWeight', 'bold');
title('Phase Portrait with Direction Field and Trajectories', 'FontSize', 18, 'FontWeight', 'bold');

% Plot equilibrium points
plot([0, 1], [0, 0], 'ro', 'MarkerSize', 12, 'MarkerFaceColor', 'red', ...
     'MarkerEdgeColor', 'white', 'LineWidth', 2);

% Add colorbar
cb = colorbar;
cb.Label.String = 'P(u,v)';
cb.Label.FontSize = 14;
cb.Label.FontWeight = 'bold';

% Set limits and grid
xlim([-1.5, 2]);
ylim([-1.5, 1.5]);
grid on;
set(gca, 'GridColor', [0.8, 0.8, 0.8], 'GridAlpha', 0.5);
set(gca, 'FontSize', 14);

hold off;

%% Method 6: Interactive Trajectory Plotting
figure('Position', [300, 300, 1200, 500]);

% Two different systems side by side
systems = {'Gradient', 'Hamiltonian'};
for sys = 1:2
    subplot(1, 2, sys);
    
    % Define vector field based on system type
    if sys == 1
        % Gradient system
        DU = -(F1*U_coarse + F2*U_coarse.^2);
        DV = V_coarse;
        title_str = 'Gradient System: \nabla P flow';
    else
        % Hamiltonian system
        DU = -V_coarse;
        DV = -(F1*U_coarse + F2*U_coarse.^2);
        title_str = 'Hamiltonian System: Symplectic flow';
    end
    
    % Contour plot
    contour(U_fine, V_fine, P_fine, contour_levels, 'LineWidth', 1.5);
    hold on;
    
    % Vector field
    magnitude = sqrt(DU.^2 + DV.^2);
    scale_factor = 0.8;
    quiver(U_coarse, V_coarse, DU./magnitude*scale_factor, DV./magnitude*scale_factor, ...
           0, 'Color', 'blue', 'LineWidth', 1.2, 'MaxHeadSize', 0.4);
    
    % Multiple trajectory starting points
    colors = {'red', 'green', 'magenta', 'cyan', 'yellow'};
    start_points = [-1.2, -0.5; -0.8, 1.0; 0.5, -1.0; 1.5, 0.8; 1.8, -0.5];
    
    for i = 1:size(start_points, 1)
        % Numerical integration of trajectories
        tspan = [0, 10];
        ic = start_points(i, :);
        
        % Define ODE system
        if sys == 1
            odefun = @(t, y) [-(F1*y(1) + F2*y(1)^2); y(2)];
        else
            odefun = @(t, y) [-y(2); -(F1*y(1) + F2*y(1)^2)];
        end
        
        [t, sol] = ode45(odefun, tspan, ic);
        
        % Plot trajectory
        plot(sol(:,1), sol(:,2), 'Color', colors{i}, 'LineWidth', 2.5);
        
        % Mark starting point
        plot(ic(1), ic(2), 'o', 'Color', colors{i}, 'MarkerSize', 8, ...
             'MarkerFaceColor', colors{i}, 'MarkerEdgeColor', 'black');
        
        % Add arrow at the end to show direction
        if length(sol) > 10
            arrow_start = sol(end-5, :);
            arrow_end = sol(end, :);
            arrow_vec = arrow_end - arrow_start;
            arrow_vec = arrow_vec / norm(arrow_vec) * 0.1;
            quiver(arrow_start(1), arrow_start(2), arrow_vec(1), arrow_vec(2), ...
                   0, 'Color', colors{i}, 'LineWidth', 3, 'MaxHeadSize', 0.8);
        end
    end
    
    % Plot equilibrium points
    plot([0, 1], [0, 0], 'ko', 'MarkerSize', 10, 'MarkerFaceColor', 'black', ...
         'MarkerEdgeColor', 'white', 'LineWidth', 2);
    
    % Formatting
    axis equal;
    grid on;
    xlabel('u', 'FontSize', 12, 'FontWeight', 'bold');
    ylabel('v', 'FontSize', 12, 'FontWeight', 'bold');
    title(title_str, 'FontSize', 14, 'FontWeight', 'bold');
    xlim([-1.5, 2]);
    ylim([-1.5, 1.5]);
    colormap(jet);
    
    hold off;
end

%% Method 7: Example usage of advanced function (defined at end of file)
% Uncomment the lines below to use the advanced function:
% create_advanced_phase_portrait(1, -1, 'gradient', 'Advanced Phase Portrait');
% create_advanced_phase_portrait(1, -1, 'hamiltonian', 'Advanced Phase Portrait');

%% Tips for Trajectory Visualization:
% 1. Use ode45 or ode23 for smooth trajectory integration
% 2. Try different starting points to reveal flow patterns
% 3. Use different colors for each trajectory
% 4. Add arrows to show direction of flow
% 5. Consider both gradient and Hamiltonian formulations
% 6. Adjust integration time to capture complete trajectories

%% ========================================================================
%% FUNCTION DEFINITIONS (Must be at the end of the file)
%% ========================================================================

function create_phase_portrait(F1, F2, plot_title)
    % Function to create a customized phase portrait
    
    % Define grid
    u = linspace(-2, 2, 800);
    v = linspace(-2, 2, 800);
    [U, V] = meshgrid(u, v);
    
    % Calculate P
    P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;
    
    % Create figure
    figure('Position', [50, 50, 600, 600]);
    
    % Create contour plot
    contour_levels = linspace(min(P(:)), max(P(:)), 15);
    [C, h] = contour(U, V, P, contour_levels, 'LineWidth', 1.5);
    
    % Styling
    colormap(jet);
    axis equal;
    grid on;
    xlabel('u', 'FontSize', 14);
    ylabel('v', 'FontSize', 14);
    title(plot_title, 'FontSize', 16, 'FontWeight', 'bold');
    
    % Add equilibrium points
    if F1 ~= 0 && F2 ~= 0
        eq_u = [0, -F1/F2];
        eq_v = [0, 0];
        hold on;
        plot(eq_u, eq_v, 'ro', 'MarkerSize', 8, 'MarkerFaceColor', 'red');
        hold off;
    end
    
    colorbar;
end

function create_advanced_phase_portrait(F1, F2, system_type, plot_title)
    % Advanced function to create phase portraits with trajectories
    % system_type: 'gradient' or 'hamiltonian'
    
    % Fine grid for contours
    u_fine = linspace(-2, 2, 800);
    v_fine = linspace(-2, 2, 800);
    [U_fine, V_fine] = meshgrid(u_fine, v_fine);
    P_fine = (1/2)*F1*U_fine.^2 + (1/3)*F2*U_fine.^3 - (1/2)*V_fine.^2;
    
    % Coarse grid for vectors
    u_coarse = linspace(-2, 2, 25);
    v_coarse = linspace(-2, 2, 20);
    [U_coarse, V_coarse] = meshgrid(u_coarse, v_coarse);
    
    % Define vector field
    if strcmp(system_type, 'gradient')
        DU = -(F1*U_coarse + F2*U_coarse.^2);
        DV = V_coarse;
    else % hamiltonian
        DU = -V_coarse;
        DV = -(F1*U_coarse + F2*U_coarse.^2);
    end
    
    % Create figure
    figure('Position', [100, 100, 900, 700]);
    
    % Filled contours
    contourf(U_fine, V_fine, P_fine, 25, 'LineStyle', 'none');
    hold on;
    
    % Contour lines
    contour(U_fine, V_fine, P_fine, 15, 'LineColor', 'black', 'LineWidth', 1);
    
    % Vector field
    magnitude = sqrt(DU.^2 + DV.^2);
    quiver(U_coarse, V_coarse, DU./magnitude*0.8, DV./magnitude*0.8, ...
           0, 'Color', 'white', 'LineWidth', 1.5, 'MaxHeadSize', 0.4);
    
    % Multiple colored trajectories
    trajectory_colors = lines(8);
    start_points = [-1.5, -1.5; -1.5, 1.5; 1.5, -1.5; 1.5, 1.5; ...
                    -0.5, 0; 0.5, 0; 0, -1; 0, 1];
    
    for i = 1:size(start_points, 1)
        % ODE system
        if strcmp(system_type, 'gradient')
            odefun = @(t, y) [-(F1*y(1) + F2*y(1)^2); y(2)];
        else
            odefun = @(t, y) [-y(2); -(F1*y(1) + F2*y(1)^2)];
        end
        
        % Solve ODE
        [t, sol] = ode45(odefun, [0, 15], start_points(i, :));
        
        % Plot trajectory
        plot(sol(:,1), sol(:,2), 'Color', trajectory_colors(i,:), 'LineWidth', 3);
        
        % Starting point marker
        plot(start_points(i,1), start_points(i,2), 'o', ...
             'Color', trajectory_colors(i,:), 'MarkerSize', 8, ...
             'MarkerFaceColor', trajectory_colors(i,:), 'MarkerEdgeColor', 'black');
    end
    
    % Equilibrium points
    if F1 ~= 0 && F2 ~= 0
        eq_points = [0, -F1/F2];
        plot(eq_points, [0, 0], 'ks', 'MarkerSize', 12, 'MarkerFaceColor', 'yellow', ...
             'MarkerEdgeColor', 'black', 'LineWidth', 2);
    end
    
    % Formatting
    colormap(parula);
    axis equal;
    xlabel('u', 'FontSize', 16, 'FontWeight', 'bold');
    ylabel('v', 'FontSize', 16, 'FontWeight', 'bold');
    title([plot_title, ' (', upper(system_type), ' System)'], 'FontSize', 18, 'FontWeight', 'bold');
    
    colorbar('FontSize', 12);
    grid on;
    xlim([-2, 2]);
    ylim([-2, 2]);
    
    hold off;
end

Answer 1

[Star Strider]

I am not certain what you want. My best guess is that you want smaller arrows in the last two figures. The fifth argument to quiver is the arrow scaling factor. I set it to 0.5 in both plots:

quiver(U_coarse, V_coarse, DU./magnitude*scale_factor, DV./magnitude*scale_factor, ...

0.5, 'Color', 'blue', 'LineWidth', 1.2, 'MaxHeadSize', 0.1);

quiver(arrow_start(1), arrow_start(2), arrow_vec(1), arrow_vec(2), ...

0.5, 'Color', colors{i}, 'LineWidth', 3, 'MaxHeadSize', 0.1);

You can of course change it to provide whatever sized arrows you want. (Setting this value to 0 turns off all scaling.) I also changed the 'MaxHeadSize' value.

Try something like this --

%% Phase Portrait Plots in MATLAB
% Creating beautiful contour plots similar to the reference images

%% Method 1: Single Plot with Your Exact Parameters
% Define the function P = (1/2)*F1*u^2 + (1/3)*F2*u^3 - (1/2)*v^2
% With F1 = 1, F2 = -1

% Define grid
u = linspace(-1.5, 2, 1000);
v = linspace(-1.5, 1.5, 1000);
[U, V] = meshgrid(u, v);

% Define parameters
F1 = 1;
F2 = -1;

% Calculate the function P
P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;

% Create figure
figure('Position', [100, 100, 800, 600]);

% Create contour plot
contour_levels = -0.4:0.1:0.4;
[C, h] = contour(U, V, P, contour_levels, 'LineWidth', 1.5);

% Enhance the plot
colormap(parula);  % You can also try: jet, hot, cool, spring, etc.
axis equal;
grid on;
grid minor;
xlabel('u', 'FontSize', 14, 'FontWeight', 'bold');
ylabel('v', 'FontSize', 14, 'FontWeight', 'bold');
title('Phase Portrait: P = \frac{1}{2}u^2 - \frac{1}{3}u^3 - \frac{1}{2}v^2', ...
      'FontSize', 16, 'FontWeight', 'bold');

% Add colorbar
colorbar('FontSize', 12);

% Set axis properties
set(gca, 'FontSize', 12, 'LineWidth', 1.2);
xlim([-1.5, 2]);
ylim([-1.5, 1.5]);

% Add equilibrium points (find where gradient = 0)
% For this system: ∂P/∂u = F1*u + F2*u^2 = u + (-1)*u^2 = u(1-u) = 0
% So u = 0 or u = 1
% ∂P/∂v = -v = 0, so v = 0
equilibrium_points_u = [0, 1];
equilibrium_points_v = [0, 0];
hold on;
plot(equilibrium_points_u, equilibrium_points_v, 'ko', 'MarkerSize', 8, ...
     'MarkerFaceColor', 'black', 'LineWidth', 2);
hold off;

Warning: Error in state of SceneNode.
String scalar or character vector must have valid interpreter syntax:
Phase Portrait: P = \frac{1}{2}u^2 - \frac{1}{3}u^3 - \frac{1}{2}v^2

%% Method 2: Multiple Subplots with Different Parameters
figure('Position', [150, 150, 1200, 400]);

% Different parameter combinations
params = [3, -3; -1, -2; 0, -2];  % [F1, F2] pairs
titles = {'(a) F_1 = 3, F_2 = -3', '(b) F_1 = -1, F_2 = -2', '(c) F_1 = 0, F_2 = -2'};

for i = 1:3
    subplot(1, 3, i);
    
    F1 = params(i, 1);
    F2 = params(i, 2);
    
    % Calculate P for current parameters
    P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;
    
    % Create contour plot
    [C, h] = contour(U, V, P, contour_levels, 'LineWidth', 1.5);
    
    % Enhance each subplot
    colormap(parula);
    axis equal;
    grid on;
    xlabel('u', 'FontSize', 12);
    ylabel('v', 'FontSize', 12);
    title(titles{i}, 'FontSize', 14, 'FontWeight', 'bold');
    
    % Set consistent axis limits
    xlim([-1.5, 2]);
    ylim([-1.5, 1.5]);
    
    % Find and plot equilibrium points for current parameters
    if F1 ~= 0
        % Solve F1*u + F2*u^2 = 0 → u(F1 + F2*u) = 0
        eq_u = [0, -F1/F2];
        eq_v = [0, 0];
        % Filter for real solutions within plot range
        valid_idx = isreal(eq_u) & eq_u >= -1.5 & eq_u <= 2;
        eq_u = eq_u(valid_idx);
        eq_v = eq_v(valid_idx);
    else
        % If F1 = 0, then F2*u^2 = 0, so u = 0
        eq_u = 0;
        eq_v = 0;
    end
    
    hold on;
    if ~isempty(eq_u)
        plot(eq_u, eq_v, 'ko', 'MarkerSize', 6, 'MarkerFaceColor', 'black');
    end
    hold off;
    
    set(gca, 'FontSize', 10);
end

%% Method 3: Enhanced Version with Better Colors and Styling
figure('Position', [200, 200, 900, 700]);

% Use a more sophisticated colormap
custom_colormap = [
    0.2, 0.1, 0.5;  % Dark blue
    0.1, 0.4, 0.8;  % Blue
    0.1, 0.8, 0.8;  % Cyan
    0.9, 0.9, 0.1;  % Yellow
    0.8, 0.4, 0.1;  % Orange
    0.5, 0.1, 0.2   % Dark red
];

% Reset to original parameters
F1 = 1;
F2 = -1;
P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;

% Create filled contour plot for better visualization
[C, h] = contourf(U, V, P, 20, 'LineStyle', 'none');
hold on;

% Add contour lines
[C2, h2] = contour(U, V, P, contour_levels, 'LineWidth', 1.2, 'LineColor', 'black');

% Apply custom colormap
colormap(custom_colormap);

% Enhanced styling
axis equal;
xlabel('u', 'FontSize', 16, 'FontWeight', 'bold');
ylabel('v', 'FontSize', 16, 'FontWeight', 'bold');
title('Enhanced Phase Portrait', 'FontSize', 18, 'FontWeight', 'bold');

% Add colorbar with custom formatting
cb = colorbar;
cb.Label.String = 'P(u,v)';
cb.Label.FontSize = 14;
cb.Label.FontWeight = 'bold';

% Plot equilibrium points
plot([0, 1], [0, 0], 'wo', 'MarkerSize', 10, 'MarkerFaceColor', 'white', ...
     'MarkerEdgeColor', 'black', 'LineWidth', 2);

% Set axis properties
set(gca, 'FontSize', 14, 'LineWidth', 1.5);
xlim([-1.5, 2]);
ylim([-1.5, 1.5]);

% Add grid
grid on;
set(gca, 'GridColor', [0.3, 0.3, 0.3], 'GridAlpha', 0.3);

hold off;

%% Method 4: Example usage of custom functions (defined at end of file)
% Uncomment the line below to use the custom function:
% create_phase_portrait(1, -1, 'Custom Phase Portrait: F_1 = 1, F_2 = -1');

%% Method 5: Phase Portrait with Direction Field and Trajectories
figure('Position', [250, 250, 1000, 800]);

% Parameters
F1 = 1;
F2 = -1;

% Define finer grid for contours
u_fine = linspace(-1.5, 2, 1000);
v_fine = linspace(-1.5, 1.5, 1000);
[U_fine, V_fine] = meshgrid(u_fine, v_fine);
P_fine = (1/2)*F1*U_fine.^2 + (1/3)*F2*U_fine.^3 - (1/2)*V_fine.^2;

% Coarser grid for vector field
u_coarse = linspace(-1.5, 2, 20);
v_coarse = linspace(-1.5, 1.5, 15);
[U_coarse, V_coarse] = meshgrid(u_coarse, v_coarse);

% Compute the vector field (gradient descent: du/dt = -∂P/∂u, dv/dt = -∂P/∂v)
% For potential systems: 
% du/dt = -∂P/∂u = -(F1*u + F2*u^2)
% dv/dt = -∂P/∂v = -(-v) = v

DU = -(F1*U_coarse + F2*U_coarse.^2);  % -∂P/∂u
DV = V_coarse;                          % -∂P/∂v

% Alternative: Hamiltonian system (symplectic)
% DU = -V_coarse;                       % ∂P/∂v  
% DV = -(F1*U_coarse + F2*U_coarse.^2); % -∂P/∂u

% Normalize vector field for better visualization
magnitude = sqrt(DU.^2 + DV.^2);
DU_norm = DU ./ (magnitude + 1e-8);
DV_norm = DV ./ (magnitude + 1e-8);

% Create filled contour plot
[C, h] = contourf(U_fine, V_fine, P_fine, 20, 'LineStyle', 'none');
hold on;

% Add contour lines
contour_levels = -0.4:0.1:0.4;
[C2, h2] = contour(U_fine, V_fine, P_fine, contour_levels, 'LineWidth', 1.2, 'LineColor', 'black');

% Add vector field
quiver(U_coarse, V_coarse, DU_norm, DV_norm, 0.5, 'Color', 'white', ...
       'LineWidth', 1.5, 'MaxHeadSize', 0.3);

% Add streamlines from various starting points
start_u = [-1.2, -0.8, -0.4, 0.4, 0.8, 1.2, 1.6];
start_v = [-1.2, -0.8, -0.4, 0, 0.4, 0.8, 1.2];
[start_U, start_V] = meshgrid(start_u, start_v);
streamline(U_coarse, V_coarse, DU, DV, start_U(:), start_V(:), [0.1, 10000]);

% Enhanced styling
colormap(parula);
axis equal;
xlabel('u', 'FontSize', 16, 'FontWeight', 'bold');
ylabel('v', 'FontSize', 16, 'FontWeight', 'bold');
title('Phase Portrait with Direction Field and Trajectories', 'FontSize', 18, 'FontWeight', 'bold');

% Plot equilibrium points
plot([0, 1], [0, 0], 'ro', 'MarkerSize', 12, 'MarkerFaceColor', 'red', ...
     'MarkerEdgeColor', 'white', 'LineWidth', 2);

% Add colorbar
cb = colorbar;
cb.Label.String = 'P(u,v)';
cb.Label.FontSize = 14;
cb.Label.FontWeight = 'bold';

% Set limits and grid
xlim([-1.5, 2]);
ylim([-1.5, 1.5]);
grid on;
set(gca, 'GridColor', [0.8, 0.8, 0.8], 'GridAlpha', 0.5);
set(gca, 'FontSize', 14);

hold off;

%% Method 6: Interactive Trajectory Plotting
figure('Position', [300, 300, 1200, 500]);

% Two different systems side by side
systems = {'Gradient', 'Hamiltonian'};
for sys = 1:2
    subplot(1, 2, sys);
    
    % Define vector field based on system type
    if sys == 1
        % Gradient system
        DU = -(F1*U_coarse + F2*U_coarse.^2);
        DV = V_coarse;
        title_str = 'Gradient System: \nabla P flow';
    else
        % Hamiltonian system
        DU = -V_coarse;
        DV = -(F1*U_coarse + F2*U_coarse.^2);
        title_str = 'Hamiltonian System: Symplectic flow';
    end
    
    % Contour plot
    contour(U_fine, V_fine, P_fine, contour_levels, 'LineWidth', 1.5);
    hold on;
    
    % Vector field
    magnitude = sqrt(DU.^2 + DV.^2);
    scale_factor = 0.8;
    quiver(U_coarse, V_coarse, DU./magnitude*scale_factor, DV./magnitude*scale_factor, ...
           0.5, 'Color', 'blue', 'LineWidth', 1.2, 'MaxHeadSize', 0.1);
    
    % Multiple trajectory starting points
    colors = {'red', 'green', 'magenta', 'cyan', 'yellow'};
    start_points = [-1.2, -0.5; -0.8, 1.0; 0.5, -1.0; 1.5, 0.8; 1.8, -0.5];
    
    for i = 1:size(start_points, 1)
        % Numerical integration of trajectories
        tspan = [0, 10];
        ic = start_points(i, :);
        
        % Define ODE system
        if sys == 1
            odefun = @(t, y) [-(F1*y(1) + F2*y(1)^2); y(2)];
        else
            odefun = @(t, y) [-y(2); -(F1*y(1) + F2*y(1)^2)];
        end
        
        [t, sol] = ode45(odefun, tspan, ic);
        
        % Plot trajectory
        plot(sol(:,1), sol(:,2), 'Color', colors{i}, 'LineWidth', 2.5);
        
        % Mark starting point
        plot(ic(1), ic(2), 'o', 'Color', colors{i}, 'MarkerSize', 8, ...
             'MarkerFaceColor', colors{i}, 'MarkerEdgeColor', 'black');
        
        % Add arrow at the end to show direction
        if length(sol) > 10
            arrow_start = sol(end-5, :);
            arrow_end = sol(end, :);
            arrow_vec = arrow_end - arrow_start;
            arrow_vec = arrow_vec / norm(arrow_vec) * 0.1;
            quiver(arrow_start(1), arrow_start(2), arrow_vec(1), arrow_vec(2), ...
                   0.5, 'Color', colors{i}, 'LineWidth', 3, 'MaxHeadSize', 0.1);
        end
    end
    
    % Plot equilibrium points
    plot([0, 1], [0, 0], 'ko', 'MarkerSize', 10, 'MarkerFaceColor', 'black', ...
         'MarkerEdgeColor', 'white', 'LineWidth', 2);
    
    % Formatting
    axis equal;
    grid on;
    xlabel('u', 'FontSize', 12, 'FontWeight', 'bold');
    ylabel('v', 'FontSize', 12, 'FontWeight', 'bold');
    title(title_str, 'FontSize', 14, 'FontWeight', 'bold');
    xlim([-1.5, 2]);
    ylim([-1.5, 1.5]);
    colormap(jet);
    
    hold off;
end

Warning: Failure at t=1.098593e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (3.552714e-15) at time t.

Warning: Failure at t=8.109230e-01. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.776357e-15) at time t.

Warning: Failure at t=2.555559e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (7.105427e-15) at time t.

Warning: Failure at t=2.237221e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (7.105427e-15) at time t.

Warning: Failure at t=6.780461e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.421085e-14) at time t.

Warning: Failure at t=4.674759e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.421085e-14) at time t.

Warning: Failure at t=5.360646e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (1.421085e-14) at time t.

●

%% Method 7: Example usage of advanced function (defined at end of file)
% Uncomment the lines below to use the advanced function:
% create_advanced_phase_portrait(1, -1, 'gradient', 'Advanced Phase Portrait');
% create_advanced_phase_portrait(1, -1, 'hamiltonian', 'Advanced Phase Portrait');

%% Tips for Trajectory Visualization:
% 1. Use ode45 or ode23 for smooth trajectory integration
% 2. Try different starting points to reveal flow patterns
% 3. Use different colors for each trajectory
% 4. Add arrows to show direction of flow
% 5. Consider both gradient and Hamiltonian formulations
% 6. Adjust integration time to capture complete trajectories

%% ========================================================================
%% FUNCTION DEFINITIONS (Must be at the end of the file)
%% ========================================================================

function create_phase_portrait(F1, F2, plot_title)
    % Function to create a customized phase portrait
    
    % Define grid
    u = linspace(-2, 2, 800);
    v = linspace(-2, 2, 800);
    [U, V] = meshgrid(u, v);
    
    % Calculate P
    P = (1/2)*F1*U.^2 + (1/3)*F2*U.^3 - (1/2)*V.^2;
    
    % Create figure
    figure('Position', [50, 50, 600, 600]);
    
    % Create contour plot
    contour_levels = linspace(min(P(:)), max(P(:)), 15);
    [C, h] = contour(U, V, P, contour_levels, 'LineWidth', 1.5);
    
    % Styling
    colormap(jet);
    axis equal;
    grid on;
    xlabel('u', 'FontSize', 14);
    ylabel('v', 'FontSize', 14);
    title(plot_title, 'FontSize', 16, 'FontWeight', 'bold');
    
    % Add equilibrium points
    if F1 ~= 0 && F2 ~= 0
        eq_u = [0, -F1/F2];
        eq_v = [0, 0];
        hold on;
        plot(eq_u, eq_v, 'ro', 'MarkerSize', 8, 'MarkerFaceColor', 'red');
        hold off;
    end
    
    colorbar;
end

function create_advanced_phase_portrait(F1, F2, system_type, plot_title)
    % Advanced function to create phase portraits with trajectories
    % system_type: 'gradient' or 'hamiltonian'
    
    % Fine grid for contours
    u_fine = linspace(-2, 2, 800);
    v_fine = linspace(-2, 2, 800);
    [U_fine, V_fine] = meshgrid(u_fine, v_fine);
    P_fine = (1/2)*F1*U_fine.^2 + (1/3)*F2*U_fine.^3 - (1/2)*V_fine.^2;
    
    % Coarse grid for vectors
    u_coarse = linspace(-2, 2, 25);
    v_coarse = linspace(-2, 2, 20);
    [U_coarse, V_coarse] = meshgrid(u_coarse, v_coarse);
    
    % Define vector field
    if strcmp(system_type, 'gradient')
        DU = -(F1*U_coarse + F2*U_coarse.^2);
        DV = V_coarse;
    else % hamiltonian
        DU = -V_coarse;
        DV = -(F1*U_coarse + F2*U_coarse.^2);
    end
    
    % Create figure
    figure('Position', [100, 100, 900, 700]);
    
    % Filled contours
    contourf(U_fine, V_fine, P_fine, 25, 'LineStyle', 'none');
    hold on;
    
    % Contour lines
    contour(U_fine, V_fine, P_fine, 15, 'LineColor', 'black', 'LineWidth', 1);
    
    % Vector field
    magnitude = sqrt(DU.^2 + DV.^2);
    quiver(U_coarse, V_coarse, DU./magnitude*0.8, DV./magnitude*0.8, ...
           0, 'Color', 'white', 'LineWidth', 1.5, 'MaxHeadSize', 0.4);
    
    % Multiple colored trajectories
    trajectory_colors = lines(8);
    start_points = [-1.5, -1.5; -1.5, 1.5; 1.5, -1.5; 1.5, 1.5; ...
                    -0.5, 0; 0.5, 0; 0, -1; 0, 1];
    
    for i = 1:size(start_points, 1)
        % ODE system
        if strcmp(system_type, 'gradient')
            odefun = @(t, y) [-(F1*y(1) + F2*y(1)^2); y(2)];
        else
            odefun = @(t, y) [-y(2); -(F1*y(1) + F2*y(1)^2)];
        end
        
        % Solve ODE
        [t, sol] = ode45(odefun, [0, 15], start_points(i, :));
        
        % Plot trajectory
        plot(sol(:,1), sol(:,2), 'Color', trajectory_colors(i,:), 'LineWidth', 3);
        
        % Starting point marker
        plot(start_points(i,1), start_points(i,2), 'o', ...
             'Color', trajectory_colors(i,:), 'MarkerSize', 8, ...
             'MarkerFaceColor', trajectory_colors(i,:), 'MarkerEdgeColor', 'black');
    end
    
    % Equilibrium points
    if F1 ~= 0 && F2 ~= 0
        eq_points = [0, -F1/F2];
        plot(eq_points, [0, 0], 'ks', 'MarkerSize', 12, 'MarkerFaceColor', 'yellow', ...
             'MarkerEdgeColor', 'black', 'LineWidth', 2);
    end
    
    % Formatting
    colormap(parula);
    axis equal;
    xlabel('u', 'FontSize', 16, 'FontWeight', 'bold');
    ylabel('v', 'FontSize', 16, 'FontWeight', 'bold');
    title([plot_title, ' (', upper(system_type), ' System)'], 'FontSize', 18, 'FontWeight', 'bold');
    
    colorbar('FontSize', 12);
    grid on;
    xlim([-2, 2]);
    ylim([-2, 2]);
    
    hold off;
end

.