clear; close all; clc; C = 1; mu = 1; rho = 1; tmax = 5; ymax = 1; ts = [0.01, 0.1]; figure(1); hold on; figure(2); hold on; for idx=1:length(ts) t = ts(idx); y = linspace(0,ymax,50); dudy = C./sqrt(mu*t/rho).*exp(-y.^2./(4*mu*t)*rho); dud…