AMATH 581 Homework 2 Shallow Liquid Simulation Erik Neumann 610 N. 65th St., Seattle, WA 98103
[email protected] November 19, 2001 Abstract A model of shallow fluid behavior is evaluated using a variety of numerical solving techniques. The model is defined by a pair of partial differential equations which have two dimensions in space and one dimension of time. The equations concern the vorticity ω and the stream function ψ which are related to the velocity field of the fluid. The equations are first discretized in time and space. The time behavior is evaluated using a Runge-Kutta ordinary differential equation solver. The spatial behavior is solved using either Fast Fourier Transform, Gaussian Elimination, LU Decomposition, or iterative solvers. The performance of these techniques is compared in regards to execution time and accuracy.
Contents 1 Introduction and Overview
2
2 Theoretical Background
4
2.1
Solving for ψ - Matrix Method . . . . . . . . . . . . . . . . . . .
5
2.2
Solving for ψ - FFT Method . . . . . . . . . . . . . . . . . . . . .
6
2.3
Discretize the Advection-Diffusion Equation . . . . . . . . . . . .
6
3 Algorithm Implementation and Development
7
3.1
Construction of Matrix A . . . . . . . . . . . . . . . . . . . . . .
8
3.2
Construction of Matrix B . . . . . . . . . . . . . . . . . . . . . .
9
3.3
Pinning the Value of ψ(1, 1) . . . . . . . . . . . . . . . . . . . . .
10
1
3.4
Comparing Solvers . . . . . . . . . . . . . . . . . . . . . . . . . .
11
3.5
An FFT problem . . . . . . . . . . . . . . . . . . . . . . . . . . .
11
4 Computational Results
11
4.1
Results for various initial conditions . . . . . . . . . . . . . . . .
11
4.2
Running times . . . . . . . . . . . . . . . . . . . . . . . . . . . .
17
4.3
Accuracy of solvers . . . . . . . . . . . . . . . . . . . . . . . . . .
17
4.4
Symmetry of Solution . . . . . . . . . . . . . . . . . . . . . . . .
18
4.5
Time Resolution Needed . . . . . . . . . . . . . . . . . . . . . . .
18
4.6
Mesh Drift Instability . . . . . . . . . . . . . . . . . . . . . . . .
18
5 Summary and Conclusions
21
A MATLAB functions used
22
B MATLAB code
23
B.1 evhump.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
23
B.2 evrhs.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
23
B.3 wh.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
23
B.4 fr.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
23
B.5 ev2.m . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
24
C Calculations
1
34
Introduction and Overview
We consider the governing equations associated with shallow fluid modeling. The intended application is the flow of the earth’s atmosphere or ocean circulation. The model assumes a 2-dimensional flow, with not much movement up or down. Another assumption is that the fluid is shallow, ie. that the vertical dimension is much smaller than the horizontal dimensions. The velocity field is given by the set of vectors ~v at each point with components u ~v = v (1) w 2
where u is the x component of the velocity, v is the y component of the velocity and so on. The