Maybe something like this:
%normal modes of string section 2.25 Stein/Wysession
% clear; close all; clc
L = 20; %string length (m)
v = 3; %string wave velocity (m/s)
xs = 10; %src x-crd (m)
dx = 0.1; %spacing of x-grid
dt = 0.1; %spacing of t-grid
tmax = 6; %lime length
tau = 0.1; %source pulse width (s)
psca = 0.5; %plot scale
x = 0:dx:L; %discrete x
nx = length(x);
t = 0:dt:tmax; %discrete time
nt = length(t);
n = 1:50; %mode integers
nn = length(n);
fprintf('nt=%d nx=%d nn=%d\n\n', nt, nx, nn)
w = n * pi * v / L; %discrete angular eigen-frequencies
F = exp( -(w*tau).^2 / 4); %source pulse
u = zeros(nt,nx,nn);
uu = zeros(nt,nx);
%%
for i = 1 : nt %time
for j = 1 : nx %space
for k = 1: nn %modes done in serially inside triple loops
u(i,j,k) = sin( n(k) *pi*xs/L)*F(k) * sin( n(k) *pi* x(j) /L)*cos( w(k) *t(i) );
end
uu(i,j) = sum( u(i,j,:) ); %sum over modes to get uu space/time matrix
end
end
%%
subplot(1,3,1)
u_sum = permute(sum(u,1),[2 3 1]);
scale = 1;
u_sum_norm = scale*u_sum./max(abs(u_sum),[],'all');
u_sum_norm_shift = u_sum_norm+(1:nn);
plot(x,u_sum_norm_shift)
set(gca,'ydir','normal','TickDir','out');
xlabel('x(m)');
ylabel('modes (n)');
title('mode')
ylim([1 nn])
%%
subplot(1,3,2)
scale = 1;
uu_norm = scale*uu./max(abs(uu),[],'all');
uu_norm_shift = uu_norm+(1:nt).';
plot(x,uu_norm_shift.','k')
set(gca,'ydir','normal','TickDir','out');
xlabel('x(m)');
ylabel('time (s)');
title('wavefield')
ylim([1 nt])
%%
subplot(1,3,3)
% set(gcf,'position', [0 50 150 200]);
colormap('turbo')
imagesc(x,t,uu);
colorbar;
set(gca,'ydir','normal','TickDir','out');
xlabel('x(m)');
ylabel('time (s)');
title('wavefield')
