source: anuga_work/development/carrier_greenspan/paul_guard/cg.asv @ 5649

Last change on this file since 5649 was 5649, checked in by duncan, 16 years ago

A start in the carrier_greenspan validation. Working on it here initially so it doesn't get released until it works.

File size: 1.6 KB
Line 
1% To run, open matlab, so to this directory and enter cg
2
3clear
4f=figure;
5slope=0.1;              %bed slope
6T=20;                   %period
7l=slope*9.8*T^2/pi^2;   %characteristice length (lo in C&G58)
8A=0.5;                  %amplitude
9dd=0.1;                 %resolution of grid in the sigma-alpha plane
10
11%define the potential phi according to the C&G solution (eqn 2.21 C&G58)
12for sigma=[1:200];
13for lambda=[1:200];
14Phi(sigma,lambda)=A*besselj(0,(sigma-1)*dd)*cos((lambda-1)*dd);
15end
16end
17
18%solve for x, t, u and eta (equations 2.16-2.19 C&G58)
19for sigma=[2:199];
20for lambda=[2:199];
21u(sigma,lambda)=-(Phi(sigma+1,lambda)-Phi(sigma-1,lambda))/dd/2./(sigma*dd);
22eta(sigma,lambda)=1/4.*(Phi(sigma,lambda+1)-Phi(sigma,lambda-1))/dd/2-1/2.*(u(sigma,lambda)).^2; % stage
23t(sigma,lambda)=1/2.*(lambda*dd)+u(sigma,lambda);
24x(sigma,lambda)=-1/4*(Phi(sigma,lambda+1)-Phi(sigma,lambda-1))/dd/2+1/16*(sigma*dd).^2+1/2*u(sigma, lambda).^2;
25end
26end
27
28%plot using Matlab interpolation routines
29figure
30ti=-50:250;
31tj=1:99;
32[XI,YI] = meshgrid(ti/10,tj/10);
33Z = griddata(x(:,:),t(:,:),eta(:,:),XI,YI);
34mesh(XI,YI,Z)
35
36%show animation using dimensional variables according to C&G58 scaling
37figure(f)
38a=0;
39for tp=1:99
40    plot(l/slope-l*XI(1,:),l-l*slope*XI(1,:),'k')
41    hold on
42    plot(l/slope-l*XI(1,:),l+l*slope*Z(tp,:),'b')
43    text(l,l+l*slope/8,num2str(tp/10*T/pi))
44    axis([0 1.1*l/slope l-l*slope/4 l+l*slope/4])
45    pause(0.1)
46    hold off
47    a=a+1;
48%     fname=[num2str(a) '.jpg'];
49%     saveas(f,fname)
50    F(a) = getframe;
51end
52movie(F,5)
53% figure
54% plot3(x, t, eta,'.')
55
Note: See TracBrowser for help on using the repository browser.