RUSSIAN JOURNAL OF EARTH SCIENCES VOL. 8, ES6004, doi:10.2205/2006ES000216, 2006
[60] The problem of sensitivity of applied tsunami problems to hydrodynamic accuracy of mathematical models and to the computational accuracy of discrete algorithms that realize the models becomes the first in the numerical solution of such problems. In the first part of this study we used the example of the problem of wave transformation in "wash-tub" model region, which was first considered by [Gusyakov and Chubarov, 1985], to analyze the peculiarities of different mathematical models reproducing the wave transformation in the regions with different depths and wave interactions with different types of boundaries (open and rigid reflecting boundaries).
[61] In the second half of the work, we study the problem of choosing the computational algorithms, which provide the required accuracy of the modeling results and formulate a number of test problems for estimating the adequacy of the developed algorithms and programs.
[62] Computational algorithms used in the problems of tsunami modeling can be divided, for example, into the algorithms designed for the solution of differential equations, algorithms of data processing (field, historical, etc.), and others. Here, we shall consider the algorithms of the first group based on finite difference methods. They are designed for calculation of wave amplitudes, velocities or transports, characteristics of the interaction with the coast, fields of maximal heights, and other functionals.
[63] The main requirements to the algorithms are directed to reach the needed accuracy, economical properties, and flexibility. Sometimes these requirements correlate well with each other but sometimes hardly solved contradictions arise. Such properties as stability margin, conservatism, invariant features, monotonicity, and lack of high frequency fluctuations facilitate gaining the satisfied accuracy in the need of using rough spatial discretization. Due to the approximation, computational algorithms inherit the parameters of the mathematical models ("hydrodynamic" and "mathematical") and introduce their own "computational" parameters, among which the most important are discretization parameters (size of computational cells, basic functions), which determine to a great degree internal properties of the algorithm, as well as additional parameters, which are used to control the scheme dispersion and dissipation, and other properties, in particular those mentioned above. We should also note the situation, which appears when asymmetric algorithms are used (multi-step splitting methods, methods of the predictor-corrector type, etc.) when special alternation of individual stages of the algorithm is required to avoid this asymmetry. In this case, symmetry is gained not at each time step.
[64] A transition to a discrete region of variation in spatial variables and time can be performed using different methods: by introducing a regular of irregular grid, which can be uniform or not uniform, either triangle, rectilinear, or curvilinear, fixed or varying in time, adjusted to the geometry of the calculation region or approximating it to a certain degree. The specific geometrical feature of calculation areas in tsunami problems, which reflects the specific features of the corresponding real basins, manifests itself in their multiply connected property and extremely complicated form of internal and external boundaries determined by the complexity of the coastline contours. In some situations, when the formulation of the problem requires the computation of the zone of inundation, such boundaries become mobile, which significantly complicates the problem. This led to the development of algorithms, which allow for the so-called continuous computation without distinguishing the water-land boundary (a type of the method of dummy regions).
[65] The efficiency of the computational algorithm depends strongly on the correctness of specifying the numerical boundary conditions. The conditions that correspond to the physics of the studied phenomenon are formulated at the moment when the mathematical model is being chosen for reproducing such processes as reflection from a wall, free pass, runup of waves on the shore including the moving one (landslide above water), or arrival of the wave into the computational area. It is desirable to approximate these "physical" boundary conditions so that the order of the accuracy of the algorithm as a whole is not changed. Usually, the specific property of the numerical method requires formulation of additional conditions, whose number and form are determined jointly by the algorithm and mathematical model. A theoretical investigation of the correctness of computational boundary conditions is possible only in the simplest cases, while in the solution of specific problems their experimental study is emphasized, which requires thoughtful tests and a large amount of computational experiments.
[66] Maintenance of the required accuracy of calculations can be gained using two methods, at least. The first of them is related to the use of more detailed input information. A quantitative measure in this case would be a spatial size, for example the size of the grid if regular discretization is used. The key problem here is the accuracy of input information, which consists of the sets of bathymetric data (digital charts), data on the roughness of the bottom and land topography, etc. A transition to a more detailed bathymetry introduces the effect of increasing accuracy, but here, one has to pay attention to the fact that more detailed charts of real basins are the result of applying certain interpolation algorithms to the initial data measured at the final number of points. Application of such procedures is inevitably accompanied by the accumulation of errors. Uniform decreasing of the grid in the entire computational area increases the resource capacity of the calculation and very frequently is not needed. An alternative to this approach is application of non-uniform grids of different type when condensation that provides increased accuracy of the results is performed, for example, in the zone of shallow depths or in the zone of their sharp variations. Such solution of the problem is accompanied by certain difficulties, for example by the fact that the computational area is multiply connected or by the behavior of the solution at the boundaries between the fragments of the computational areas of different size of the grid.
[67] Application of the algorithms of increased order of approximation can be considered as the second method of increasing the accuracy. As a rule, the choice of computational methods of such class complicates significantly the algorithm, increases the number of points involved in the calculation, and complicates the formulation and realization of the boundary conditions. This aspect of providing the accuracy of calculations requires thorough grounds and detailed study to make a reasonable solution.
[68] Increasing of the accuracy of numerical modeling can be gained also by more accurate specification of the boundaries in the calculation area. This approach becomes most pressing in the case, when the objective of modeling is determination of the characteristics of wave regime at the boundary and near boundary points. Different methods exist to gain such increase in the accuracy, for example, by application of curvilinear grids adjusted to the contour of the boundary or by transition from traditional finite difference schemes to the methods of finite element type, finite volumes, and others, which use triangle or irregular quadrangular grids. Certain difficulties for the finite element methods in tsunami problems are related to the reproduction of a number of characteristic boundary conditions.
[69] The problem of providing the computational stability in tsunami problems has its own specific features. First of all, one should differentiate the local and global problems in time. As a rule, stability of calculations can be gained by standard methods in the solutions of local problems designed for modeling short-term processes occurring in limited basins (first of all by using theoretically stable methods that have reliable stability margins and mechanisms, which filter out non-physical high frequency oscillations). Long-term modeling, in the course of which the wave can propagate over transoceanic distances, requires additional software tools, permanent monitoring of the state of calculated fields with the objective of early distinguishing of potential sources of instability and their timely removal either automatically or as a result of interactive interruption of computational expert. The development of instability during long calculations is usually related to the so-called "nonlinear instability". The time of its occurrence and spatial localization can be hardly performed in advance or this is practically impossible. In the test problems considered by the authors we made an attempt to correlate these effects with certain peculiarities of bottom topography. However no universal criterion has been developed yet.
[70] Stability of "short term" calculations depends strongly on the smoothness of bottom topography in the area considered in the study and smoothness of the initial data, which are input data. They are determined, for example, using other computational algorithms or as a result of processing of field measurements. The computational experience of the authors indicates that such dependence manifests itself to a certain degree when mathematical models are used, whose equations include high order derivatives (non-linear dispersion models, account for elastic effects at the free surface in modeling the interaction of waves with floating platforms, airdromes, etc.). The problem of smoothness of the input information can be solved using special smoothing algorithms, which to a certain degree of sophistication include the mechanisms of artificial dissipation. If necessary, these algorithms can be used to remove high-frequency computational oscillations of the grid scale. Such work requires high computational qualification. If the work is performed by a dilettante, it frequently leads to the attempts to interpret purely computational effects as physical ones. An example of such type of errors would be given below. We note, however that the determination whether one or another effect in the results of numerical modeling is "physical" requires that each problem should be solved using sequentially decreasing grids. Only conservation of the found effect at any scale of discretization allows us to start interpreting it as an existing phenomenon. Naturally, we have to mention the correct choice of the time step and a possibility of its variation from one layer to another in the solution of nonlinear problems with attempts to maintain the "Courant number" as close to unity as possible. The solution of the problem of computational stability (prevention of the collapse of numerical solution) in the case of long-term calculation requires special thoroughness in the formulation of the modeling problem. For example, the mechanisms of smoothing and artificial dissipation, or others mentioned above could be an admissible means of providing stability if it is needed to reproduce the long-term dynamics of the most general characteristics of the phenomenon and determine general regularities in the variation of wave parameters in space and time. A different case is in the situations when it is required to determine absolute values of wave heights, velocities, and runup values. Such calculations should be preceded by thorough methodological computational experiments over real topography using the grids of different size, various computational algorithms, and different stabilization mechanisms.
[71] The considerations described above would be illustrated on the examples related to the solution of a number of test problems including those considered in the previous sections about wave transformations in the "wash-tub" basin.
Citation: 2006), Principles of numerical modeling applied to the tsunami problem, Russ. J. Earth Sci., 8, ES6004, doi:10.2205/2006ES000216.
Copyright 2006 by the Russian Journal of Earth Sciences (Powered by TeXWeb (Win32, v.2.0).