How to 'crop' a graph?

Question refers to the last section of code, which produces a graph of U vs time. I want to know how to remove the first cycle of that graph. Thanks in advance!

%% Solves 2 PDE's for ETA and U using finite differencing 
    % (centralised with lagged friction)

close all;
clear all;              % Initial Clearing 
clc;

% Discretise Domain
n = 21;                 % Number of eta points in our domain

disp(['Number of elevation points in the grid = ',num2str(n)])

Channel_Length = 22000; % m (thus 20km in length) we will have three 
                        %  time-levels per space discritisation because of 
                        %  the time filter employed (which requires t-1:t+1 
                        %  for each i point.
% Set up values
eta(1:n,1:2) = 0;
eta(1:n,3) = NaN;
u(1:n-1,1:2) = 0;       % Velocity staggard grid (eta drives velocity)
u(1:n-1,3) = NaN;

% Set up domain
d = zeros(1,n)+10;      % Create coinstant bathyemtry at 10m 
b(1:n) = 180;           % Channel width (note spatially constant)

% Set grid size and timestep 
dx = Channel_Length/(n-1);    % Length/n of point = dx

% Domains
xeta = 0:dx:(n-1)*dx;         % x-values for plotting eta
xu = dx/2:dx:(n-1)*dx-dx/2;   % x-values for plotting u

%
dt = 20;    % Timestep
c = dt/dx;

% CFL check for time-step stability
if dt >= dx/sqrt(9.81*min(d))
    disp('Warning!: dt is too large')
end
tmax = 8000;            % Number of timesteps in simulation
time = 0:dt:tmax*dt;    % Sim_time in seconds
cd = 0.001;             % Bottom drag coefficient. **Change this!**

% To get the hydrodynamics stated in the assignment
disp(['The bottom drag coefficent is = ',num2str(cd)])
r = 0.05;               % Constant for Asselin filter 

%% Initialise variables and set initial conditions

h1 = 1.9;       % Amplitude of boundary at point 1
g1 = 180;       % Phase of boundary at point 1
hn = 2.6;       % Amplitude of boundary at point n
Phase_Diff = ((24/60)/12.4224)*360;
gn  =  g1+Phase_Diff;    % Phase of boundary at point n 
                        %  12 degrees more than other boundary
% H1 = 0.0125*5; G1 = 100;
% Hn = 0.0750*5; Gn = 050;
g = 9.81;
count = 1;  % Counter for output of eta and u history
nans = 50;  % Output frequency (animation and u/eta history)
            %  Hint: (nans*dt*20)/3600
Eta = nan(1+tmax/nans,1);   % Precondition output arrays
U = nan(1+tmax/nans,1);     %           "

% Store intial conditons
Eta(count) = eta(1,1);
U(count) = 0.5*(u(5,1)+u(6,1));     % Store initial condition
count = count+1;

%% Computation with Central Finite Difference and Lagged Scheme

for it = 1:tmax                 % Computational loop
    
    % Time setting
    t = it*dt;
    tim = t/3600;
    tc1 = tim/12.4224;
    tc2 = tim/6.2112;
    sim_tim(count) = tim;       % For plotting
    
    % Elevation boundary conditions
    eta(1,3) = h1*cos(2*pi*tc1-g1*pi/180);  % +H1*cos(2*pi*tc2-G1*pi/180);
    eta(n,3) = hn*cos(2*pi*tc1-gn*pi/180);  % +Hn*cos(2*pi*tc2-Gn*pi/180);
    
    % Update elevations
    for i = 2:n-1
        eta(i,3) = eta(i,1)-c*((b(i)*(d(i)+eta(i,2))+b(i+1)*...
            (d(i+1)+eta(i+1,2)))*u(i,2)-(b(i-1)*(d(i-1)+eta(i-1,2))...
            +b(i)*(d(i)+eta(i,2)))*u(i-1,2))/b(i);
    end
    
    % Update velocities
    for i = 1:n-1
        db1 = -g*(eta(i+1,2)-eta(i,2))/dx;
        db2 = -cd*(u(i,1)*abs(u(i,1)))/((0.5*(d(i)+eta(i,1)+d(i+1)+...
            eta(i+1,1)))^(4/3));
        u(i,3) = u(i,1)+2*dt*(db1+db2);
    end
    
    % Time smoothing (Asselin filter)
    for i = 1:n
        eta(i,2) = eta(i,2)+0.5*r*(eta(i,1)-2*eta(i,2)+eta(i,3));
    end
    
    for i = 1:n-1
        u(i,2) = u(i,2)+0.5*r*(u(i,1)-2*u(i,2)+u(i,3));
    end

    % Store results and plot (animation)
    if rem(it,nans) == 0
        subplot(2,1,1)
            plot(xeta,eta(:,3),'-ob','MarkerSize',5,'LineWidth',1.5);
            title(num2str(tim,'%4.1f'))
            ylabel('\eta (m)');
            ylim([-2.5 2.5]);xlim([0 xeta(n)]);
            subplot(2,1,2)
        plot(xu,u(:,3),'-ob','MarkerSize',5,'LineWidth',1.5);
            ylabel('u (m s^-^1)');  xlabel('distance (m)')
            ylim([-2 2]);xlim([0 xeta(n)]);
        drawnow
        % Try to output the velcoity at the middle of the channel (note: Is
        % this correct?)
        U(count) = mean([u(10,3),u(11,3)]);     % Output velocity history
        E(count) = mean([eta(10,3),eta(11,3)]); %           "
        Ust(count) = u(1,3);
        eta1(count) = eta(1,3);
        Ue(count) = u(end-1,3);
        eta22(count) = eta(end-1,3);
    	u(:,count) = u(:,3); % Why 3?  =  =  for time-filtering 
        count = count+1;
    end
    
    % Shuffle (update time indices)
    eta(:,1) = eta(:,2); eta(:,2) = eta(:,3);
    u(:,1) = u(:,2); u(:,2) = u(:,3);

end

% Plot u time history
figure
    Tim = 0:dt*nans:tmax*dt;
    plot(Tim,U,'-ob','MarkerSize',5,'LineWidth',1.5);
        hold on
    plot([0 max(Tim)],[0 0],'--b','MarkerSize',5,'LineWidth',1.5);
        title('U'); 
        xlabel('Time (s)'); 
        ylabel('U Mid Channel (m/s)');

Answer 1

[Star Strider]

I’m not certain what you want.

One solution:

% Plot u time history
zci = @(v) find(v(:).*circshift(v(:), [-1 0]) <= 0);                    % Returns Approximate Zero-Crossing Indices Of Argument Vector
Uzxidx = zci(U);
figure
    Tim = 0:dt*nans:tmax*dt;
    plot(Tim(Uzxidx(3):end),U(Uzxidx(3):end),'-ob','MarkerSize',5,'LineWidth',1.5);
        hold on
    plot([0 max(Tim)],[0 0],'--b','MarkerSize',5,'LineWidth',1.5);
        title('U'); 
        xlabel('Time (s)'); 
        ylabel('U Mid Channel (m/s)');

This finds the zero-crossings and omits everything up to the third (at ‘Uzxidz = 46’).

Another option is to use the xlim function:

        xlim([Tim(46) max(xlim)])

Both of these work, producing slightly different results.

[Dan Williams]

Legend, thanks! It was as simple as xlim, couldn't get it to work :/

[Star Strider]

As always, my pleasure!

I didn’t see xlim in your code, so I can’t determine the problem you were having wih it.