Czech-Japanese Seminar in Applied Mathematics 2006


DNAPL Flow In Heterogeneous Porous Media - A Multi-Scale Stochastic Approach

K. Barnhart1, T. H. Illangasekare1, D. W. Dean2 and T. F. Russell2
1Center for Experimental Study of Subsurface Environmental Processes (CESEP), Colorado School of Mines, Golden, Colorado
2Department of Mathematics, University of Colorado at Denver, Denver, Colorado


This research explores the use of stochastic differential equations (SDEs) to model multiphase flow in heterogeneous aquifers, specifically the flow of DNAPLs in satu-rated soils. A fundamental assumption used in the model formulation is that the flow of fluid particles is described by a stochastic process and that the positions of the fluid particles over time are governed by the law of the process. Darcys law and the continuity equation are used to derive a Fokker-Planck equation from these expressions. The Ito calculus is then applied to derive a SDE for the non-wetting phase. Stochastic models like these are typically used to model dispersion processes. Such models, in their usual form, cannot represent multiphase barrier effects which occur at the interface between adjacent materials of different permeabilities. In this model, the control of the flow of DNAPL particles across the interface is accomplished using a jump term which derives from the Ito formula. The jump term is based on capillary diffusivity and the pressure-saturation curves of the sands forming the interface. An additional microscale addition is presented which is intended to simulate irregulari-ties in the pore structure at the interface. From empirical data, lognormal distributions of grain size are used to generate a probabilistic percolation at the continuum scale. These same probability distributions may also be used to simulate fingering due to microheterogeneities in a (supposedly) homogeneous media.

Numerical Simulation of Air Flow and Pollution Transport in the Atmospheric Boundary Layer

P. Bauer

Czech Technical University in Prague, Czech Republic.


We create a mathematical model based on Navier-Stokes equations for media flow and diffusion-convection equation describing pollution transport, and solve the model using finite element methods (FEMs). We use a simple algebraic model with turbulent viscosity. The FEMs allows us to use different terrain shapes and examine the terrain influence. Essential part of the work is numerical analysis, which examines the properties of algorithms with respect to numerical parameters, and analyzes the selected cases of air flow in the boundary layer of atmosphere in the case of instant, steady or periodically interrupted source of pollution.

Cardio MRI Data Segmentation using PDE of Allen-Cahn Type

R. Chabiniok

Charles University, Prague and Institute for Clinical and Experimental Medicine, Prague, Czech Republic.


The work deals with segmentation of image data using the algorithm based on numerical solution of the geometrical evolution partial differential equation of the Allen-Cahn type. This equation has origin in the description of motion by mean curvature. The equation has diffusive character. The diffusion process can be controlled by the segmented signal, so that the edges of the objects or areas can be found. The method is applied to the problem of automatic segmentation of the left heart ventricle from the cardio MRI images. The segmentation is a necessary step for the estimation of significant indicators of the myocardial function such as ejection fraction or kinetic parameters of the wall thickening during cardiac cycle - parameters that describe clinical situation of the myocardium.

Keywords: image processing, segmentation, PDE, Allen-Cahn equation, cardio MRI

Existence, uniqueness and L2(Ω) a priori estimates for discrete solution of nonlinear tensor diffusion in image processing

O. Drblíková

Slovak University of Technology in Bratislava, Faculty of Civil Engineering, Department of Mathematics and Descriptive Geometry, Slovakia.


We introduce two-dimensional semi-implicit diamond cell finite volume scheme for solving the nonlinear tensor diffusion in image processing. We focus on the proof of existence and uniqueness of the discrete solution given by our scheme. The main idea of this proof is based on a bounding the gradient in tangential direction by the gradient in normal direction. At the end, L2(Ω) a priori estimates for discrete solution of our scheme are proved.

Finite Element Method applied to the Geodetic Boundary Value Problem

Z. Fašková

Slovak University of Technology in Bratislava, Faculty of Civil Engineering, Department of Mathematics and Descriptive Geometry, Slovakia.


The poster deals with the external geodetic boundary value problem for the determination of the height anomaly. We apply Finite Element Method to this problem with Neumann and Dirichlet boundary conditions and compare the solution with the quasigeiodal heights generated directly from EGM96 (the Earth gravitational model 1996) based on spherical harmonics with completely different mathematical strategy.

Flux-based level set method on unstructured grids

P. Frolkovič

Simulation Technology Center, University of Heidelberg, Germany.


New method for numerical solution of nonlinear degenerate level set equation will be presented. The method is formally second order accurate and consistent when applied with 2D/3D unstructured grids using different types of elements. It was implemented using a multilevel grid structure ("multi-grids") by applying a local refinement and coarsening of the computational mesh in time. The method is closely related to standard "flux-based" finite volume methods and can be implemented by relatively minor changes to existing codes for nonlinear advection-diffusion problems. Several 2D/3D examples will be shown involving ones with topological changes of the interface and with a grid adaptivity. Some applications will be discussed including the modeling of the growth of an atherosclerotic lesion or a transport model coupled with a flow of two immiscible fluids.

Benchmark Solutions for the Two-Phase Porous-Media Flow

R. Fučík

Czech Technical University in Prague, Czech Republic.


The contribution presents a discussion of one- and two-dimensional models of the two-phase flow through porous media used in verification of more complex numerical models. We give information how a particular solution of the simplified model based on the McWorther-Sunada formula can be generalized. To enlarge the class of admissible boundary and initial conditions, we offer a numerical algorithm solving the transport equation for saturation, which is based on the Finite-Difference Method in space and time and yields values of the solution at given time levels and on a spatial grid of positions. The use of the algorithm is demonstrated on a series of computations in one-dimensional spatial domain.

Application of Multiple-Precision Arithmetic to Numerically Unstable Problems

H. Fujiwara

Graduate School of Informatics, Kyoto University, Japan.


We make direct approach in numerical computation of some numerically unstable problems using multiple-precision arithmetic. Discretization of an ill-posed problem is numerically unstable in most cases. Numerical instability leads the rapid growth of computational errors then computation fails.

Some kinds of computational errors are treated in numerical analysis and we focus on the rounding error in the research. Rounding error comes from discretization of the real number and its arithmetic on digital computers, and we approximate them with about 15 digits in the standard computations. We propose the use of multiple-precision arithmetic which enables us to approximate real numbers with desired accuracy to control rounding errors.

A multiple-precision arithmetic software, which is designed and implemented by the authors for fast and large scale scientific computations, is also introduced. Proposing software works with the programming language C++ or Fortran90 on personal computers, workstations, and a supercomputer and mainly implemented in their assembly languages.

Morphological instability of a pure material in a mixed phase

P. Guba

Comenius University Bratislava, Faculty of Mathematics, Physics and Informatics, Slovakia.


The stability of a mixed-phase, or mushy, region formed during the directional solidification of a pure substance is investigated. The model adopted incorporates a non-equilibrium phase-change transition and diffusive transport of heat in the absence of convection. Two prototype problems are considered: the growth of a mixed phase into a thermally non-undercooled liquid, and that into an undercooled one. The dimensionless governing equations are parametrized in terms of a Stefan number λ and a phase-relaxation number τ. Asymptotic solutions for the disturbance growth rate are derived in the limits of short and long wavelengths. It is shown that the mixed-phase system is linearly stable in the case of a non-undercooled liquid, and is unstable in the case of an undercooled liquid, at all wavelengths and irrespective of λ and τ. Comparison with the corresponding results for the classical Stefan problem clarifies the status of the mushy-region regularisation.

Comparison of Several Finite Difference Methods for Magnetohydrodynamics in 1D and 2D

P. Havlík, R. Liska

Czech Technical University in Prague, Faculty for Nuclear Science and Physical Engineering, Czech Republic.


Magnetohydrodynamics (MHD) is a system of conservation laws obtained as an extension of Euler equations for ideal gas dynamics by including magnetic pressure and equations for magnetic field evolution. The structure of the system and its solution is richer than that of Euler equations and consists in 1D of seven Riemann waves including Alfven, fast and slow magnetoacoustic waves. We compare composite, central, WENO and CWENO finite difference schemes and numerical methods from packages Flash and Nirvana. In 1D, we present several Riemann problems, intermediate shock formation problem and smooth periodic problem which served as order of accuracy test problem as its analytical solution exists. For 2D test problems we have chosen tests regularly appearing in the literature. They include the fast rotor problem, the Orszag-Tang vortex problem, the blast problem and the shock cloud interaction problem. In general, it seems that the best results are obtained by Flash and WENO, although WENO was unable to compute the blast and the shock cloud interaction problem. Finally, we are interested in generalization the numerical methods for MHD system into cylindrical r-z geometry.

The software used in this work was in part developed by the DOE-supported ASC / Alliance Center for Astrophysical Thermonuclear Flashes at the University of Chicago

Computer Assisted Analysis for the existence of home/heteroclinic orbits

Y. Hiraoka

Hiroshima University, Department of Mathematical and Life Sciences, Japan.


In this talk, I would like to present a new approach to prove the existence of homo/heteroclinic orbits in vector field. The approach uses exponential dichotomy properties to rigorously calculate Melnikov integrals.

Why the Classical Nucleation Theory (CNT) Cannot Work

R.F. Holub1, P.K. Hopke2 and M. Beneš3
1Colorado School of Mines, Department of Environmental Science and Engineering, Department of Physics, Golden, USA.
2Clarkson University, Department of Chemistry, Potsdam, USA.
3Czech Technical University in Prague, Faculty for Nuclear Science and Physical Engineering, Czech Republic.


There is a general consensus that nucleation - and bubble formation - are plagued by variuos problems. One problem is so called capillary approximation. The well known bulk surface tension and free energies that are used to explain the formation of the critical embryo do not really work. The problem appears to be that, while macroscopically correct, the dynamics of the process, which has to go through a "nanoscopic", critical phenomena phase, behaves fundamentally different. A review article "Recent developments in the kinetic theory of nucleation" [1], uses single molecule master equation, however, even this approach does not appear to explain why they obtain somewhat higher nucleation rates than CNT. There appears to be a basic "unquantifiability" of this phenomenon, similar to the surface diffusion conundrum, as we have previously reported [2]. We are suggesting, therefore, that by invoking a reproducible - but to some extent random - anomalous 'geoaerosol' phenomenon we call "'very fast' transport of particles over large distances...", and by assuming it is of a quantum nature (a kind of "3-D surface diffusion", see slide 7 in [3]), we can potentially overcome the CNT problems.

[1] E. Ruckenstein and Y.S.Djikajev, "Advances in Colloid & Interface Science", 118, 51-72, (Dec 30, 2005)

[2] R.F.Holub, M.Beneš and B.D.Honeyman, "Interfacial phenomena in fluid dynamics: linking atomistic and macroscopic properties: can they explain the transport anomalies?", Proceedings of Czech Japanese Seminar in Applied Mathematics 2004, pp 57-62,

[3] R.F.Holub, lectures given at CVUT. FJFI, June 2006,

Numerical Simulation of Inductively Coupled Plasma Flows Using the Object-Oriented Solver COOLFluiD

R. Honzátko1, H. Deconinck2, D. Vanden Abeele2, A. Lani2, M. Panesi2 and T. Quintino2
1Czech Technical University in Prague.
2von Kármán Institute for Fluid Dynamics, Rhode-St-Genese, Belgium.


When a shuttle or capsule enters a planetary atmosphere from space, it faces very high temperatures, reaching several thousands degrees $K$. To prevent damage to the spacecraft, it needs to be protected by a heat shield, made from thermal protection materials. Samples of these new materials may be tested in high-temperature wind tunnels in which they are subjected to a hot jet of inductively heated thermal plasma.

Since high pressure inductively coupled plasmas are characterized by low Mach numbers, they should be considered as incompressible in the sense that density hardly varies due to pressure variations in the flow field. Nevertheless, density still changes due to the strong temperature variations and chemical reactions in the plasma. Therefore, a strictly incompressible and low Mach number flow solver in finite-volume method is considered. It is validated against several benchmark solutions.

The effect of an electromagnetic field generated by coils surrounding a plasma torch is incorporated.

The work deals with flows under the assumption of a local thermodynamic equilibrium with a fixed elemental composition, e.g. 21% of nitrogen and 79% of oxygen.

Numerical results of inductively coupled plasma flows under operational conditions of a plasmatron used at von Kármán Institute are presented and compared to the solutions computed by the former code developed at von Kármán Institute. Both a qualitative and quantitative comparison is made.

Keywords: inductively coupled plasma, local thermodynamic equilibrium, finite-volume method

Challenges in the Modeling the Behavior of Multiphase Flow in Naturally Heterogeneous Porous Media

T. Illangasekare

Colorado School of Mines, Department of Environmental Science and Engineering, Golden, USA.


Remediation of soil and water contaminated with organic chemicals that are in the form of non-aqueous phase liquids is a significant and a challenging problem in the environmental field. The physical problems of relevance involve the flow and entrapment behavior of multiphase fluids and the mass transfer that occurs from the non-aqueous to the flowing aqueous phase. Numerical models that capture these basic processes are used to obtain a fundamental understanding and in the design and evaluation of strategies for cleanup. This paper discusses some of the challenges that are faced in the development and validation of these models in applications dealing with complex natural systems with geologic heterogeneities.

On solutions to a BVP related to Sobolev-Poincare inequalities

T. Idogawa

Department of Electronic and Information Systems, Faculty of Systems Engineering, Shibaura Institute of Technology, Japan.


We consider a BVP (|ux|p-2ux)x+|u|q-2u=0 in (a,b) with u(a)=u(b)=0. This problem is related to the Sobolev-Poincare type inequality. We will talk about the regularity of solutions and the behavior of them w.r.t. parameters p and q.

On a reaction-diffusion model for smolder patterns under micro-gravity

K. Ikeda

Mathematical Institute, Tohoku University, Japan.


Various finger-like smoldering patterns under micro-gravity were observed in experiments. For theoretical understanding of such pattern phenomena, a model of reaction-diffusion system has been proposed. In this talk, we prove the existence and uniqueness of a solution for this reaction-diffusion system. We also consider a large-time behavior of solutions.

Finite difference schemes for Landau-Lifshitz equation

T. Ishiwata

Department of Mathematics, Faculty of Education, Gifu University, Japan.


We are concerned with finite difference schemes for Landau-Lifshitz equation which describes the evolution of spin fields in nonequilibrium continuum ferromagnets.

This equation has length-preserving property and energy conservation (Heisenberg case) or dissipation property. It is well-known that the conventional finite difference scheme does not inherit these properties. In this talk, we will propose some schemes which inherit these properties and show some numerical results.

This talk is based on a joint work with A. Fuwa and M. Tsutsumi.

Turbulence in the environmental flow

Z. Jaňour

Institute of Thermomechanics AS CR, Czech Republic.


Processes in the atmospheric boundary layer are complicated owing to combined effects of the Earth?s rotation, buoyancy forces, surface drug forces, and complex surface topography. It is extremely difficult to solve this synergetic problem analytically and experiments in situ are expensive and give only particular results (especially in complicated topography). The main source of information for solving of this problem is simulation. Therefore the system of the equations of motions for geophysical flow is introduced; turbulence modelling in the atmosphere is discussed and some examples of the simulation are demonstrated.

Computer simulations of laser beam interaction with foams

T. Kapin, R. Liska

Czech Technical University in Prague, Faculty for Nuclear Science and Physical Engineering, Czech Republic.


Low-density foam layers may be utilized in inertial fusion targets for smoothing of transverse inhomogeneities in ablation pressure. Here we present results of one- and twodimensional simulations of laser interaction with solid foams considering both approaches of foam modelling, ie., structured and uniform models. Computer implementation is based on Lagrangian and Arbitrary Eulerian - Lagrangian (ALE) approach for one- and twodimensional studies, respectively. According to theoretical predictions, ablative and explosive regimes of laser propagation through foam that are dependent on mean density and pore thickness of foam can be observed in simuations. Qualitative comparison with recent experimental results is performed.

Some gauge invariant estimates for Ginzburg-Landau equations

H. Kasai

Departmant of Mathematics and Infomatical Science, Fukushima University, Japan.


The Ginzburg-Landau (GL) equations is a model equation for superconductivity.This equation has various types of solutions called Meissner solution, n-vortex solution and trivial solution. Since, GL equations has invariance under infinite dimensional transformation called gauge transformation, it is not trivial to define stable state without some constraint condition (gauge fixing).

In this talk, we will consider some apriori estimates without gauge fixing. For this sake, some equations of quantity which have gauge invariance are considered. These apriori estimates can give information to the dynamics of vortices or an interface between normal state and superconducting state.

Exponential decay of the first eigenvalue of an elliptic problem with large drift

M. Kimura

Faculty of Mathematics, Kyushu University, Japan.


This study is concerned with the asymptotic behaviour of the first eiginvalue of an elliptic problem with a large drifting term. The magnitude of the drift is controled by a parameter. Under a suitable condition, it is known that the first eigenvalue tends to zero exponentially if the parameter tends to infinity. We show some numerical examples and prove some basic behaviours of the eigenvalue in multidimensional case. We also prove a precise behaviour of such exponential smalless of the eigenvalue in one dimensional case.

Extensions of the Mizukami-Hughes method for convection-diffusion equations

P. Knobloch

Charles University, Faculty of Mathematics and Physics, Prague, Czech Republic.


The Mizukami-Hughes method was proposed in [1] for the numerical solution of scalar two-dimensional steady convection-diffusion equations using conforming triangular piecewise linear finite elements. Several improvements of this method were recently introduced in [2]. A distinguished feature of the Mizukami-Hughes method is that it satisfies the discrete maximum principle and gives very accurate discrete solutions in convection-dominated regime. In the talk we shall discuss various extensions, generalizations and modifications of this method. In particular, we shall concentrate on the application of the method to convection-diffusion-reaction equations, three-dimensional problems and other types of finite elements.

[1] A. Mizukami, T.J.R. Hughes, A Petrov-Galerkin finite element method for convection-dominated flows: An accurate upwinding technique for satisfying the maximum principle, Comput. Methods Appl. Mech. Eng. 50 (1985) 181-193.

[2] P. Knobloch, Improvements of the Mizukami-Hughes method for convection-diffusion equations, Preprint No. MATH-knm-2005/6, Faculty of Mathematics and Physics, Charles University, Prague, 2005.

Numerical solution of unsteady flow using artificial compressibility method

K. Kozel1, P. Louda1, J. Příhoda2

1Dept. of Technical mathematics, FME, Czech Technical University.
2Institute of Thermomechanics, Czech Republic.


The work presents an implicit artificial compressibility method applied to Navier-Stokes equations for steady as well as for unsteady flows. Two modifications of unsteady numerical solution by finite volume method are considered. First one uses large artificial compressibility parameter and the iterative solution approximates unsteady evolution of flow. The second approach introduces dual time (artificial time) into the system of equations, where first time serves as iterative time only and the second one is physical time. Later approach is found more robust and reliable. The computational results are shown for laminar flow over circular cylinder and turbulent synthetic jet flows.

Numerical solution of higher and lower Mach number flows

K. Kozel, P. Punčochářová, J. Fürst, J. Horáček

Dep. of Technical Mathematics, Faculty of Mechanical Engineering, Czech Technical University in Prague, Czech Republic.


The work deals with numerical solution of two compressible fluid problems. Firstly, authors considered steady flows in DCA 8% cascade for increasing upstream Mach numbers M∈(0.813,1.13). The cascade flows were suggested in Institute of Thermomechanics (IT CAS CZ) by Mr. Dvorak and flows were investigated experimentally. The structure of flow seems to be very complicated. There is possible to observe subsonic and supersonic part, shock wave structure, interaction of shock wave and boundary layer etc. We investigated these flows numerically using finite volume method for system of Euler equations. These numerical results are compared to experimental results of IT CAS CZ using comparison of several regimes with increasing upstream Mach numbers.

The second problem that is considered in this work is unsteady flows in a 2D channel where one part of the channel wall is changinhg as a given function of time. Authors considered two types of the channel: an unsymmetrical and a symetrical channel. The flow has very low upstream Mach numbers M∈(0.02,0.1) and is described by the system of Navier-Stokes equations for compressible (laminar) flows. The problem is numerically solved by MacCormack finite volume scheme. Moved grid of quadrilateral cells is considered in the form of conservation laws using Arbitrary Lagrangian-Eulerian method.

Non-Equilibrium Thermodynamics Conception of Energy and Mass Transfer in Living Tissues and in Fuel Cells

F. Maršík, S. Převorovská, V. Klika, O. Mičan

Institute of Thermomechanics, Academy of Sciences of the Czech Republic.


Non-equilibrium thermodynamics some times called as irreversible thermodynamics is capable to explain not only the processes for many nonliving systems, but also to discover some new phenomena and relations for living tissues and for biological systems at all. In the connection with the numerical simulation this theory is the effective tool for the more intense understanding of the processes running in complex systems.

The muscle performance is demonstrated on the hemodynamical model of human cardiovascular system, which allows to evaluate many hidden parameters, e.g. mechanical powers of heart muscles, the influence of dialysis and supporting pumps on heart loads, etc. Such simulation results in more accurate diagnostics of the patient's disorders.

As a second example is studied the remodellation of the bone structure under mechanical stimuli. The ability of bone to adapt its shape and its internal composition gives some insight into the complexity of this phenomenon. It is shown that one consequence of dynamical loading is the formation of cortical bone from cancellous bone.

Electrochemical processes running in fuel cells are challenge for phenomenological description, especially, when cross effects like diffusion of oxygen and water are affected by electroosmosis and capillarity. Competition between the concentration diffusion of oxygen and the electroosmosis results in the anode dehydration. Including the water transport by the capillarity can this defect suppress. The appropriate transport of oxygen and water can be improved by some combination of number hydrophilic pores with a corresponding number of small hydrophobic pores.

Computational methods in image analysis

K. Mikula

Slovak University of Technology in Bratislava, Faculty of Civil Engineering, Department of Mathematics and Descriptive Geometry, Slovakia.


In many applications computers analyse images or image sequences quality of which can be poor, e.g., they are contaminated by a noise and/or boundaries of image objects are partly missing (e.g. in medical imaging, in scene with occlusions or ilusory contours). We will discuss how nonlinear partial differential equations can be used to denoise and segment such images, complete the missing boundaries, etc. Dicretizations by the variational methods of the geometrical image selective smoothing and segmentation equations (Riemannian mean curvature flow and Perona-Malik problem) will be presented. The computational results in bio-medical image processing will be given.

Modelling Multiphase Flow in Heterogeneous Porous Media

J. Mikyška

Czech Technical University in Prague, Czech Republic.


Multi-phase models that simulate the behavior of non-aqueous phase liquids in porous media can be used to obtain fundamental understanding of the complex behavior and predict the fate of waste chemicals in the subsurface. Existing models have limitations in simulating highly heterogeneous systems to be able to represent realistic field conditions. The presentation reports development of a new multiphase flow code called VODA. It starts with a brief introduction of the mathematical model of the multiphase flow in porous media. Then, the Control Volume Finite Element (CVFE) discretization is described and finally, several numerical experiments are presented concerning validation of the developped code. Experimental convergence analysis is carried out for two well-known one- dimensional two-phase flow problems. Examples of several further two-phase flow computations in a heterogeneous medium will be also given.

Numerical simulation of dislocation dynamics

V. Minárik, J. Kratochvíl, K. Mikula, M. Beneš

Czech Technical University in Prague, Czech Republic, Slovak Technical University, Bratislava.


The aim of this contribution is to present the current state of our research in the field of numerical simulation of dislocations moving in crystalline materials. The simulation is based on recent theory treating interactions of dislocation curves and dipolar loops, both occurring in the material and interacting by means of forces of elastic nature. The mathematical model describes the motion and interaction laws for one dislocation curve and finite number of dipolar loops placed in 3D space. The interactions occur not only between the dislocation curve and dipolar loops but also between dipolar loops themselves, which makes the model more complex. Equations of motion for a parametrically described dislocation curve are discretized by the flowing finite volume method in space. The interaction force is computed for each dipolar loop and along the discretized curve. The resulting system of ordinary differential equations is solved by a higher-order time solver.

Relaxation oscillation of seven point vortices

T. Nakaki

Graduate School of Sciences, Hiroshima University, Japan.


The motion of assembly of point vortices in a plane is one of classical and interesting problems in the field of fluid mechanics. Periodic, quasi-periodic and chaotic motions are known, for example. Relaxation oscillation, which is a special case of (quasi-)periodic motions, is one of our subjects.

This oscillation consists of steady and rapid stages, which occur alternatively. In mathematical point of view, such an oscillation is induced by heteroclinic orbits. The existence of heteroclinic orbits is already studied for n=3,4,5, where n is the number of point vortices. In this talk, we show that heteroclinic orbits exist when n=7. Some numerical simulations are demonstrated. We also touch on the cases where n=9,11,....

Numerical verification of solutions for heat convection problems

M. T. Nakao

Faculty of Mathematics, Kyushu University, Japan.


We present several results on cmputer assisted approaches for solutions of the two-dimensional Rayleigh-Benard convection problems. First, we will describe on a basic concept of our numerical verification method to prove the exsistence of the steady-state solutions based on the infinite dimensional fixed-point theorem using Newton-like operator with the spectral approximation and the constructive error estimates. Next, we show some verification examples of several exact non-trivial solutions for the given Prandtl and Rayleigh numbers. Furthermore, a computer assisted proof of the existence for a symmetry breaking bifurcation point will be presented. Furthermore, we will also mention about the extension of these results to the three dimensional problems with some prototype numerical results.

Extended Fast Marching Method for the Sail Distance

T. Nishida

Department of Mathematical Informatics, the University of Tokyo, Japan.


Suppose that a flow field is given in a two or three dimensional space, and consider a vehicle sailing in the flow. If there is no flow of water, the vehicle can move in any direction at the same maximum speed. If the water flows, on the other hand, the speed of the vehicle is anisotropic; the vehicle can move faster in the same direction as the flow, while it move only slowly in the direction opposite to the flow direction. Thus, a sail distance is defined as the smallest time necessary for the vehicle to move from one point to another.

This distance appears in technological and utilitarian situations. An airplane in a windy sky is one of the examples. When flights from Tokyo to each city of European countries are planned, wind such as westerly is the most important in the factors which have effects on the flight route, and hence we should take the sail distance into the consideration. In addition, ships on ocean current encounter the same problem.

It is easy to define the sail distance, but it is difficult to compute the sail distance. If we know the path achieving the smallest time, we can calculate the distance. However, we cannot know it in advance, and hence it is hard to calculate the distance directly.

In order to overcome this difficulty, we construct a numerical method. For this purpose, we reduce the problem of computing the sail distance to a boundary value problem of a partial differential equation. This idea is the same as the idea for reducing the problem of computing the Euclidean distance to a boundary value problem of the eikonal equation. When we set the speed of the vehicle to 0, then our proposed equation is identical with eikonal equation, and hence our equation is the generalization of eikonal equation.

For the eikonal equation, the fast marching method is well-known as an efficient and stable computational technique, and we would like to apply this method to our equation. However, the fast marching method does not work well.

In this talk, we propose a new numerical method by extending the scheme of the fast marching method and show the efficiency and the stableness of the proposal method by numerical experiments.

On a numerical scheme for the Willmore flow

T. Oberhuber

Czech Technical University in Prague, Czech Republic.


The contribution deals with the Willmore flow of graphs [1] which is described by the relaxation of a functional depending on the mean curvature. We describe a numerical scheme based on the method of lines. The spatial derivatives are discretized by the finite-difference method on a uniform grid and the resulting system of ODEs is solved by the Runge-Kutta method. We deliver several qualitative examples of computation using the derived scheme as well as numerical studies of convergence demonstrating ability of the scheme.


1. K. Deckelnick, G. Dziuk, "Error estimates for the Willmore flow of graphs", Preprint Nr. 28/2004, Universität Magdeburg, 26 pp. (2004), to appear in Interfaces and Free Boundaries

Bifurcations for Turing's instability without SO(2) symmetry

T. Okuda

Division of Mathematical Science, Graduate School of Engineering Science, Osaka university, Japan.


Bifurcation structures of Turing's reaction diffusion system with mixed boundary conditions are studied.

We don't a priori know the eigenfunctions for the linearized problem since mixed boundary conditions break the SO(2) symmetry.

We will discuss how the neutral stability curves change and , as a result, how the bifurcation diagrams change by modifying boundary conditions from usual Neumann to Mixed one .

Droplet motion with some contact angle

S. Omata

Kanazawa University, Division of Mathematical and Physical Sciences, Graduate School of Natural Science & Technology, Japan.


We investigate motion of droplet on a surface. Driving force comes from contact angle and surface tensions. The problem become hyperbolic or parabolic free boundary problem with volume conservation constraint. We have derived equations and developed a numerical method to solve this problem. An approximate solutions are also constructed.

Computer analysis of fractal sets

P. Pauš

Czech Technical University in Prague, Czech Republic.


The work discusses the results obtained by the approximate computer numerical methods applied in the analysis of dimension for selected class of fractal sets. We present several modifications of the box-counting method, show successful results as well as difficulties encountered during their development. We also deliver several theoretical results concerning the theoretical aspects of the topic.

On the dispersion error of finite element meshes

J. Plešek, R. Kolman, M. Okrouhlík

Institute of Thermomechanics CAS, Czech Republic.


Spatial discretization of a continuum introduces dispersion error to a numerical solution of wave propagation. In the introduction, review is made of fundamental approaches used to derive the truncation error of the finite element method in a dynamic analysis, valid for elements with linear shape functions.

In the second part, recent results accomplished by the authors are summarized, namely the extension of dispersion theory to quadratic finite elements, following the lines of reasoning introduced by Belytschko and Mullen [1] for one-dimensional elements and those of Abboud and Pinsky [2], concerning the scalar Helmholtz equation. The main conclusions drawn may be stated as follows. i) A spurious optical branch in the spectrum existed. ii) The associated modes possessed infinite phase velocity, finite group velocity and strongly focused polarization. iii) It was further shown, in terms of dispersion curves, that the quadratic elements had much more favourable properties than the linear ones. This may, however, not be true of diagonalized (lumped) mass matrices.

The paper closes with the discussion of procedures suitable for mass lumping, which is essential for explicit integration methods to work efficiently. While the mass matrices of linear elements can be lumped rather easily by merely summing up the row coefficients, the same method does not work properly with higher order elements since it creates negative diagonal members. There were some remedies proposed, namely the Hinton-Rock-Zienkiewicz (HRZ) lumping scheme. As an illustrative example, the dispersion induced by the HRZ diagonal matrix was investigated. It turned out that the HRZ method spoiled the otherwise good dynamic properties of quadratic finite elements. Other approaches had better be chosen.

Acknowledgments to GACR 101/06/0914 and GAAV A100100626 under AV0Z20760514.

1. T. Belytschko, R. Mullen, "On dispersive properties of finite element solutions," In: Modern Problems in Elastic Wave Propagation, Wiley, New York 1978.
2. N. N. Abboud, P. M. Pinsky, "Finite element dispersion analysis for the three-dimensional second-order scalar wave equation," International Journal for Numerical Methods in Engineering, v. 35, p. 1183--1218, 1992.

Numerical Solution of Flows in Atmospheric Boundary Layer

K. Seinerová, K. Kozel, L. Beneš, T. Bodnár

Czech Technical University in Prague, Czech Republic.


The work deals with a numerical solution of viscous flows in the atmospheric boundary layer. The mathematical model is based on the system of Navier-Stokes equations. Numerical solution uses the finite difference method. Numerical results over 2D and 3D single hill are presented and discussed.

Finite Element Analysis of 200 Million DOF Pressure Vessel model with ADVENTURE system on the Earth Simulator

R. Shioya

Faculty of Engineering, Kyushu University, Japan.


We have been developing an advanced general-purpose computational mechanics system, named ADVENTURE, which is designed to be able to analyze a three dimensional finite element model of arbitrary shape over 100 million Degrees Of Freedom (DOF) mesh. The one of main process modules for solid analysis, named ADVENTURE_Solid, is based the hierarchical domain decomposition parallel algorithm and employs the balancing domain decomposition as a solution technique for linearized equations. The ADVENTURE_Solid has been successfully implemented on a single PC, PC clusters and MPPs with high parallel performances.

In order to realize a virtual demonstration test, the solid module is implemented on the Earth Simulator with minor modification for a vector processor and applied for an implicit dynamic elastic analysis such as seismic response of a precise nuclear pressure vessel model of 35 million DOF mesh with 128 nodes (1,024 Processing Elements), and succeeded in solving 300 unsteady steps in about 4.3 hours.

Furthermore, the system is applied for one time step of implicit dynamic elastic analysis of a nuclear pressure vessel model of 200 million DOF mesh about 13 minutes on the 128 nodes (1,024 PEs), and archived 2.1 TFLOPS, which is 25.8% of the peak performance.

Keywords: parallel finite element method, balancing domain decomposition, seismic response analysis, ADVENTURE system

Numerical simulation of reaction-diffusion dynamics

R. Straka

Czech Technical University in Prague, Czech Republic.


Reaction--diffusion system Brusselator is theoretical model of nonlinear chemical reaction in stirred reactor tank. In the reaction scheme the initial components A and B are transformed into products D and E via the reaction intermediates X and Y. It is well known that this model exhibit various types of solution depending on the characteristic length of reactor tank L. Keeping L small L∈<0,0.513) we obtain trivial invariant set - stable fixed point. If L increases we obtain stable symmetic periodic solution - we crossed first Hopf bifurcation point L*=0.513. The branch of stable antisymmetric solutions bifurcate in neighbourhood of L+=1.225 and further develop to invariant torus. Development of this torus were studied using methods of Poincare mapping and its fractal dimension and Lyapunov characteristic exponents. Symmetry of solutions and period of periodic solution were studied as well.

Volume-Preserving Motion

K. Švadlenka

Kanazawa University, Division of Mathematical and Physical Sciences, Graduate School of Natural Science & Technology, Japan.


We study parabolic and hyperbolic PDE with volume-conservation constraint originating in more complex free-boundary problems. The general overview of the research is presented focusing on the volume constraint and its numerical computation.

3-dimensional Thermal Convection Computations in the Glass Production Furnace

D. Tagami

Dep. of Intelligent Machinery and Systems, Faculty of Engineering, Kyushu University, Japan.


A class of finite element schemes for nonstationary thermal convection problems with temperature-dependent coefficients are considered. Because of an argument based on the energy method, optimal error estimates for the velocity and the temperature without any stability conditions have been established. A model of the practical melting glass convection is then computed by a finite element method based on such a scheme mathematically justified. Some numerical results are demonstrated to investigate differences of the energy efficiency among some data configurations.

Cardio MRI - dynamic examination of the heart hemodynamic functions

J. Tintěra and R. Chabiniok

Institute of Clinical and Experimental Medicine, Prague, Czech Republic.


The presentation summarizes current potential of magnetic resonance imaging (MRI) to visualize and quantify hemodymic function of the heart as the ejection fraction, stroke volume, myocardial wall flexibility and morphology. Dynamic cardio MRI provide also data for the modeling of the heart motion, myocardial wall thickening and blood flow pattern in normal or pathological cases.

Three general MRI methods as data sources are presented: a) dynamic CINE imaging, b) CINE with radio-frequency (RF) tagging grid and c) flow velocity mapping with phase-contrast principle. First method is widely used to segment heart ventricle volumes in time and also to segment myocardial wall to assess its thickness during cardiac cycle as one of the most important parameters of the heart function. For the proper segmentation (in best way as automatic procedure) we require high image contrast between blood pool and myocardial tissue, sufficient spatial resolution and also good time resolution (the fastest techniques can already provide images with less then 10 ms per time frame). Some techniques to achieve such high quality MRI data will be shown. Second method - RF tagging - is mostly used to measure stress and strain within myocardium. Some models were already developed to visualize strain lines in the left ventricle myocardial wall from these MRI data. Third presented method is used to quantify blood flow velocity to monitor flow pattern in heart ventricles and valves or in large connecting vessels (aorta).

An objective of this survey presentation is to motivate theoretical work in modeling of the liquid and tissue dynamics and possibly correlate these models in their capability to predict pathological evolution in bio-material with real "in vivo" imaging of the human tissue.

A numerical method for Gauss curvature flow based on crystalline algorithm

T. Ushijima

Department of Mathematics, Faculty of Science and Technology, Tokyo University of Science, Japan.


The Gauss curvature flow in R3 makes a smooth strictly convex surface shrink with the outward normal velocity equals to the Gauss curvature with negative sign. For this problem, we constructed an approximation algorithm using so-called crystalline motion of polyhedrons [3].

The crystalline motion was introduced by Taylor [2] and Angenent & Gurtin [1] to analyze crystal growth mathematically. The most typical crystalline motion of polygon in R2 makes each edge of a polygon keep the same direction but move with the normal speed inversely proportional to its length. It is known that such motion can be used to approximate the curvature flow of plane curves. In [3], we extended such motion of polygon to a motion of polyhedron and proved that it approximated the Gauss curvature flow.

In this talk, we will consider a numerical method for the Gauss curvature flow based on this approximation algorithm.

This is a joint work with H. Yagisita.


[1] S. Angenent and M. E. Gurtin, Multiphase thermomechanics with interfacial structure 2. Evolution of an isothermal interface, Arch. Rational Mech. Anal., 108 (1989), 323--391.
[2] J. E. Taylor, Constructions and conjectures in crystalline nondifferential geometry, Differential Geometry} (eds. B. Lawson and K. Tanenblat), Pitman, (1991), 321--336.
[3] T. K. Ushijima and H. Yagisita, Convergence of a three-dimensional crystalline motion to Gauss curvature flow, Japan Journal of Industrial and Applied Mathematics, 22 (2005), 443--459.

Asymptotic behavior of solutions to an area-preserving motion by crystalline curvature

S. Yazaki

Faculty of Engineering, University of Miyazaki, Japan.


Asymptotic behavior of solutions of an area-preserving crystalline curvature flow equation is investigated. In this equation, the area enclosed by the solution polygonal curve is preserved, while its total interfacial energy keeps on decreasing. In the case where the initial curve is an admissible convex polygon, from Brünn and Minkowski's inequality and the theory of dynamical systems, we can show that the asymptotic shape of a solution polygon is the boundary of the Wulff shape. In the case where the initial curve is an essentially admissible convex polygon or is not necessarily convex, a few results are known. The talk will be focused on the known results with sketch of the proof and ongoing research with some numerical simulations.

Numerical studies of a Two-Scale Model for Liqud-Phase Epitaxy with convection

V. Chalupecky1, Ch. Eck2 and H. Emmerich3
1Faculty of Nuclear Science and Physical Engineering, Czech Technical University in Prague, Czech Republic
2Institute for Applied Mathematics, University of Erlangen, Germany
3Institute of Minerals Engineering, RWTH Aachen, Germany


The aim of this poster is to present a two-scale model for liquid phase epitaxy [1,2] and its numerical scheme together with some numerical experiment. We consider a diffusion-convection model for transport of macroscopic concentration of atoms and temperature in a three-dimensional volume above the epitaxial surface. The microscopic growth of the surface is simulated by a Burton-Cabrera-Frank type model, which is coupled to the macroscopic quantities via coupling conditions. The numerical scheme is based on the method of lines using finite-difference discretization in space. The resulting system of ordinary differential equations is solved by a Runge-Kutta method of fourth-order. The two-scale nature of the model allows us to use a coarse grid in the three-dimensional volume and a fine grid on the two-dimensional epitaxial surface, thus making the algorithm more efficient. We also present results of some numerical experiments demonstrating the physical relevance of the model presented.

[1] Ch. Eck, H. Emmerich, Models for liquid phase epitaxy, DFG Preprint 146, 2004
[2] V. Chalupecky, Ch. Eck, H. Emmerich, Computation of nonlinear multiscale coupling effects in liquid phase epitaxy, DFG Preprint 196, 2006