- 117 Views
- Updated On :
- Presentation posted in: General

Nonclassified activity of CML RFNC-VNIITF. Oleg V. Diyankov Head of Computational MHD Laboratory. Russian Federal Nuclear Center-All-Russian Institute for Technical Physics. Snezhinsk, Chelyabinsk region (the Urals), Russia

Nonclassified activity of CML RFNC-VNIITF

Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author.While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server.

- - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - -

Nonclassified activity of CML RFNC-VNIITF

Oleg V. Diyankov

Head of Computational MHD Laboratory

Russian Federal Nuclear Center-All-Russian Institute for Technical Physics

- Snezhinsk, Chelyabinsk region (the Urals), Russia
- About 9000 employees, among which there are scientists: physicists, chemists, mathematicians; designers, engineers, etc.
- RFNC has many divisions. We’re representing the division of theoretical physics and applied mathematics.

Division of Theoretical Physics and applied mathematics

- 260 scientists in computational mathematics and theoretical physics.
- Computational MHD lab was created on the 1st of April 1996.
- 14 scientists are working in it.

- 2.5D MHD code
- 2D Irregular Grid for Mathematical Modeling
- Linear Solvers for Flows in Porous Media Modeling
- Difference Schemes for Hyperbolic Systems Treatment
- 3D Elastic-Plastic Modeling of Processes of Ceramics Formation
- 3D Gas Dynamic Code for Instability Investigation
- Development of Special Software Tools for CERN

- The Physical Model Realized in the Code
- The Application of the Code to Laser Beam Interaction with the Matter Modeling
- The Application of the Code to Liner Magnetic Compression Modeling
- The Application of the Code to the Plasma Channel Formation Modeling

2½D MHD MAG Code

O.V.Diyankov, I.V.Glazyrin,

S.V.Koshelev, I.V.Krasnogorov, A.N.Slesareva, O.G.Kotova

Russian Federal Nuclear Center - VNIITF

P.O.Box 245, Snezhinsk, Chelyabinsk Region, Russia

The main goal of MAG code creation was the necessity

of modeling of hot dense plasmas in magnetic field.

The system of equations used in MAG code is determined by Braginskii model for one-temperature case:

r

Let us assign indices x and y in the case of the axial symmetry to the r- and z- components of vector variables and r- and z- components of independent spatial variables correspondently,and index z to the - component of vector variables.

Then we receive a basic system of equations, which depends on the parameter of symmetry ( = 0 - the plane symmetry, = 1 - the symmetry is axial).

The equations have been written in arbitrary moving coordinate system, and then splitting into two systems has been performed. The diffusive terms have been splitted to a separate system of equations, the remained terms produce a quasilinear hyperbolic system.

Two equations for x- and y- components of magnetic field was written in form of an z-component vector potential A:

The first system is a hyperbolic one and it describes the ideal MHD flows in arbitrary moving coordinate system.

The second one is a diffusive system of equations. It includes the equations for energy, z (or ) – component of magnetic field and z (or ) component of vector potential.

Boundary conditions:

• applied pressure: P|b= f(t), where P|b means the pressure at the corresponding boundary, f(t) is a given time dependent function;

• rigid wall: un|b = 0, where un is a normal to boundary component of mass velocity;

• piston: un|b = f(t);

Boundary conditions:

They are: given temperature, given heat flux, heat flux as a function of temperature. These conditions may be written in the form:

Here, T|b is the temperature at the corresponding boundary, , are numerical parameters (=0 and =1 for given temperature and =1 and =0 for heat flux), f(T|b, t) is a given function. nT|b is a normal to the boundary component of the temperature gradient.

Boundary conditions:

• symmetry: Bz /n=0, B/n=0, Bn =0 ;

• conducting wall: Bz=0, B=0, Bn =0 ;

• axis: Bz=0, B=0, Bn =0, where Bis used for the axial symmetry and Bz for the plane one;

• given current:

Bz=2p j /c – plane symmetry,

B= 2p I /cr– axial symmetry,

where Bn is the normal to the boundary component of the magnetic field, j - current density, I is the whole current in z - direction, r - upper radius.

- Lagrange (no mesh reconstruction)
- Euler (grid nodes are returned to original positions at the n-th time step)
- Local (only eight neighbor nodes are used for new node coordinates determination)
- Algebraic (new nodes coordinates are calculated by bilinear interpolation of boundaries nodes coordinates in mathematic coordinate system)
- Poison equation solution is used for new coordinates determination

Equation of continuity:

Euler equations for velocity components:

x component:

y component:

z component:

Equation for A:

Equation for z component of magnetic field:

Equation for total energy:

Equation for the square root of the metric tensor:

Here uk is a projection of mass velocity vector to k-th vector of a local basis, uk=uk, where k is a covariant local basis vector, vk,Bk are determined in the same way, v - coordinate system velocity, g is a determinant of metric tensor with covariant component gij, k - mathematical coordinates.

r = const.- equation of continuity,

g = const.- equation for determinant of metric tensor,

u = const.- Euler equation for velocity.

Equation for A:

Equation for Bz:

Equation for energy:

Where:

Difference sheme for diffusive equations system:

Experimental Setup

Simulation

If plasma liner implosion is used as plasma radiation source one needs to receive the uniform plasma column. The ideal configuration for such type of radiation source could be an annular pinch (see Fig.1, left). The plasma, imploded from large radius, reaches the axis and the efficient radiator is formed. But this scheme is unsuitable because of MHD instabilities which present the great danger to uniformity of the liner implosion. So the compression achieved in corresponding experiments is substantially lower than one predicted by one-dimensional (1D) MHD calculations.

Magnetic field,

-component

[10 MGs]

Density [g/cc]

The equation of radiation transfer along a ray has the following form:

Temperature evolution,

time = 0.1, 0.3, 0.5, 0.8 ns.

Energy of laser pulse is

4 kJ/cm2. Triangle pulse with duration time of 1 ns. Focal spot 50 mm.

Temperature evolution.

time = 0.1, 0.3, 0.5, 0.8 ns.

Energy of laser pulse

4 kJ/cm2. Triangle pulse shape with duration time of 1 ns. Focal spot 50 mm.

Background plasma density was 10-6 g/cm3

The plasma is heated up to high

temperature (temperature reaches 450 eV).

A motion of the SW should be spherical but the SW dynamics

propagating along the laser pulse

direction differs from the perpendicular SW one.

As the velocity of SWs is approximately equal to 107 cm/sec, the perpendicular SW leaves the region of the laser absorption (the focal radius) after 0.3-0.5 ns.

- Gas Dynamics
- Poisson Equation Solver
- Maxwell Equation Solver
- Heat Transfer

IGM Code

Oleg V. Diyankov

Sergei S. Kotegov

Vladislav Yu. Pravilnikov

Yuri Yu. Kuznetsov

Aleksey A. Nadolskiy

RFNC – ARITPh

supported by LLNL grant B329117

- The IGM code was created to perform 2D flows modeling. The main feature of the code is the possibility of large deformations accounting.
- 3 physical processes are taken into account now: gas dynamics, heat conduction and Poisson equation.
- The main advantage before well-known finite element codes is the possibility of arbitrary deformations description.

- The flow region is initially covered by a set of Voronoi cells, and then at each time step the grid is reconstructed.
- This allows neighbor points to move free in any direction, so they may move very far from each other.
- The GUI interface for the IGM code (it is called CELLS) allows to put in initial data (geometry, matters, initial distribution of the values), and to look through the received results.

- Plasma physics (instability study, laser produced plasma, and so on).
- High velocity impact.
- Heat transfer.
- Electrostatic fields.

Benefits:

- Local orthogonality
- Local uniformity

and

difference scheme:

boundary condition:

Integrate over Voronoi cell

where

where

Poisson problemsf(r)= 10x8(x-2)8[36(x-1)2+2x(x-2)]y10(2-y)10+ 10y8(y-2)8[36(y-1)2+2y(y-2)]x10(2-x)10(x,y)=x10(2-x)10y10(2-y)10

- The source code, written in C++, has approximately 100000 lines.
- Typical calculation takes from 10 minutes up to 1 hour on the SGI R10000.
- RAM needed: 800 bytes per cell.

Nonstationary Preconditioned Iterative Methods for Large Linear Systems Solving

Oleg V. Diyankov

Vladislav Y. Pravilnikov

Natalia N. Kuznetsova

Sergei V. Koshelev

Igor V. Krasnogorov

Dmitriy V. Gorshkov

Aleksey A. Nadolsky

- ·Introduction
- ·Iterative Methods
- ·Preconditioners
- ·The PLS Code

- The work has been performed according to the contract with ExxonMobil Upstream Technology Corporation.
- The main goal was to select methods which are the most appropriate ones to the ExxonMobil specific problems.
- The results of many scientists and especially prof. Yousef Saad, Dr. Edmond Chow, prof. Kolotilina, Dr. Eremin have been used.

The goal of the work is to solve the linear system of equations

Ax=b,(1)

whereA is a large sparse matrix, bisa right hand side vector, xisa vector of unknowns.Preconditioned Iterative Methodscould be used for solving iteratively such systems.

IterativeMethods

Iterative methods can be written in the general form as follows:

(2)

where

- the sequence of nonsingular

- the sequence of real parameters.

matrices,

Or, if we denote

, formula (2) can

be expresses in the following form:

(3)

Iterative methods are called stationary, if

don’t depend upon the iteration count j. Otherwise,

iterative methods are called nonstationary.

The most known stationary iterative methods are the Jacobi method, the Gauss-Seidel method and Successive Overrelaxation (SOR) method.

The most known nonstationary iterative methods are the Conjugate Gradient (CG) method and Minimal Residual (MinRes) method and their modifications.

As a rule, for nonstationary iterative methods the solution xjon J-th iteration is searched through minimization of a quadraticfunctional from xj.

The Minimal Residual (MinRes) method

The MinRes iterations can be written in the

following form:

where

Modifications: The General Minimal Residual (GMRES) method and Quasi-Minimal Residual (QMR) method.

In general the minimal residual method can be expressed in the following form:

where Hjis a matrix, which depends from some parameters and which is chosen to minimize |b-Axj|.

This method can be applied to systems with nonsymmetrical matrix A.

The Conjugate Gradient (CG) method

The solution xj on J-th iteration is constructed

as an element of

is minimized.

so that

The convergence is guaranteed for systems with symmetric positive definite matrix A.

Modifications: The BiConjugate Gradient (BiCG)method, The Conjugate Gradient Squared (CGS) method, The BiConjugate Gradient Stabilized (BiCGStab) method.

The BiConjugate Gradient (BiCG) method

The BiCG method generate two CG-type

bi-orthogonal sequences of vectors, one based

on A, and one on AT.

This method can be used for systems with

nonsymmetrical and nonsingular matrix A.

Preconditioners

(4)

The preconditioner is called incomplete, if it is constructed on a subset of indices:

- There is three basic type of preconditioners:
- The Incomplete Cholesky (IC) factorization;
- The Incomplete LU (ILU) factorization;
- The Approximate Inverse (AI) factorization.

The Incomplete Cholesky (IC) factorization

The Incomplete LU (ILU) factorization

- The algorithm of the ILU(0) (ILU without extension of stencil) construction can be described as follows.
- For each row of A we should fulfill the following
- steps:
- For the I-th row we select all nonzeros having indices less than I;
- For each J we calculate aij = aij / ajj;
- For only nonzeros in I-th row with indices K greater than I we fulfill the eliminationaik = aik – aij * ajk.

- There is two main strategies for construction of
- a preconditioners with extended stencil:
- structural;
- fill-in with drop tolerance.
- The algorithm of the FILU (fill-in ILU with
- drop tolerance eps ) can be described as follows.
- For each row of A we should fulfill the following
- steps:
- For the I-th row we select all nonzeros having indices less than I;
- For each such J we calculate aij = aij / ajj;
- For all nonzeros in I-th and J-th rows with indices K greater than I we fulfill the elimination aik = aik – aij * ajk if aik exists and for all other nonzero ajkwe check the drop condition | aij * ajk| > eps * |Ai|.

The Approximate Inverse (AI) factorization

The algorithm of the AI construction can be

described as follows.

- For each row of matrix A we should fulfill following steps:
- The determination of the fill stencil;
- The calculation of the Mi coefficients from auxiliary linear system P * mi = ai, where miand ai consists from nonzeros of rows Miand Ai respectively.

The PLS code

Brief description of the PLS code

- The PLS (Parallel Linear Solver) is the C++ object-oriented code for solving of a large sparse linear systems. Linear systems are considered as multi-block systems. For solving of these systems we use preconditioned iterative methods.
- General steps of solution process are:
- Transformation of a multi-block system to a single block system;
- Generation of a Preconditioner;
- Execution of an Iterative Method;
- Inverse transformation of obtained solution to a multi-block structure.
- The general object-oriented program languages concepts such as the polymorphism and the derivation has been used in the code.

- We’ve realized four version of the code:
- Serial Unix version;
- Parallel Unix version based on MPI;
- Threads Unix version;
- Serial Windows-NT version.
- The PLS code can be used both as command-line utility and as library with program interface. The user can use these interfaces in his application to choose the iterative method, the preconditioner and call solver.
- We’ve implemented the following iterative methods:CG, BiCG, CGS, BiCGStab, MMR(l) (the Modified Minimal Residual method),BiCGStab(L).
- We’ve implemented the following preconditioners: ILU, MILU, FILU, MFILU, AIP, AIP(L).
- Moreover, the special solver of the block Gauss-Siedel type has been realized for systems with a 3x3 block matrix.

- The results, presented here, have been obtained on Pentium PC cluster system (PII/400).
- Equations marked as ExxonMobil ones are samples of real equations, appearing in ExxonMobil applications. Their structure and specific characteristics is ExxonMobil proprietary information.

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

0.65

BiCGStab

16.81

17.80

1

AIP(1)

7.74

BiCGStab

39.87

47.95

3

ParAIP(1)

3.14

ParBiCGStab

21.35

24.49

7

ParAIP(1)

1.35

ParBiCGStab

10.18

11.53

Task 1.

The solution of Poisson equation with discretization on irregular grid constructed on Voronoi cells.

Matrix size: N=50249, NNZ=349945.

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

9.79

BiCGStab

26.59

37.19

1

AIP(1)

17.81

BiCGStab

69.85

88.48

3

ParAIP(1)

7.57

ParBiCGStab

43.89

51.46

7

ParAIP(1)

3.13

ParBiCGStab

21.65

24.78

Task 2

ExxonMobil’s Mixed Cube problem (mc30).

Matrix size: N=110700, NNZ=593100

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

26.56

BiCGStab

102.95

131.41

1

AIP(1)

42.50

BiCGStab

237.43

281.88

3

ParAIP(1)

18.27

parBiCGStab

152.93

171.20

7

ParAIP(1)

7.79

parBiCGStab

100.67

108.46

15

ParAIP(1)

3.61

parBiCGStab

55.62

59.23

Task 3

ExxonMobil’s Mixed Cube problem. (mc40)

Matrix size: N=260800, NNZ=1406400

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

58.86

BiCGStab

242.32

304.97

1

AIP(1)

83.78

BiCGStab

722.97

810.52

3

ParAIP(1)

36.46

parBiCGStab

411.47

447.93

7

ParAIP(1)

15.71

parBiCGStab

233.10

248.81

15

ParAIP(1)

7.11

parBiCGStab

124.04

131.15

Task 4

ExxonMobil’s Mixed Cube problem. (mc50)

Matrix size: N=507500, NNZ=2747500

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

2.47

BiCGStab

31.99

35.61

1

AIP(0)

2.8

BiCGStab

271.87

275.82

1

AIP(1)

20.80

BiCGStab

114.61

136.53

1

AIP(2)

261.08

BiCGStab

144.24

406.48

3

ParAIP(1)

8.91

parBiCGStab

71.10

80.01

7

ParAIP(0)

0.62

parBiCGStab

116.81

117.42

7

ParAIP(1)

3.62

parBiCGStab

34.23

37.85

7

ParAIP(2)

50.53

parBiCGStab

49.43

99.96

15

ParAIP(1)

1.75

parBiCGStab

21.69

23.44

Task 5

ExxonMobil’s Mfem problem. (mfem_32_64_16)

Matrix size: N=134656, NNZ=704000

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

6.41

BiCGStab

139.16

149.33

1

AIP(1)

72.21

BiCGStab

810.89

887.03

3

ParAIP(1)

30.89

parBiCGStab

412.42

443.31

7

ParAIP(1)

13.23

parBiCGStab

196.56

209.79

15

ParAIP(1)

6.17

parBiCGStab

153.08

159.25

Task 6

ExxonMobil’s Mfem problem. (mfem_48_96_24)

Matrix size: N=450432, NNZ=2395008

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

18.19

BiCGStab

105.48

124.89

1

FILU(0.02, 0.2)

6.21

Geo(BiCGStab)

87.71

94.35

3

ParAIP(1)

2.84

parGeo

parBiCGStab

36.97

39.81

7

ParAIP(1)

1.40

parGeo

parBiCGStab

20.23

21.63

15

ParAIP(1)

0.74

parGeo

parBiCGStab

11.87

12.61

Task 7

ExxonMobil’s Geo problem. (geo20)

Matrix size: N=60983, NNZ=1430843

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.01, 0.1)

100.23

BiCGStab

130.88

234.18

1

FILU(0.01, 0.1)

31.76

Geo(BiCGStab)

341.21

374.44

3

ParAIP(1)

10.47

parGeo(parBiCGStab)

224.51

234.98

7

ParAIP(1)

4.80

parGeo(parBiCGStab)

105.95

110.75

15

ParAIP(1)

2.71

parGeo(parBiCGStab)

63.03

65.74

Task 8

ExxonMobil’s Geo problem. (geo30)

Matrix size: N=200073, NNZ=5203713

Nproc

Preconditioner

Tprec

Iterative Method

Tim

Ttotal

1

FILU(0.02, 0.2)

208.22

BiCGStab

125.75

337.89

1

FILU(0.02,

0.2

13.74

Geo(BiCGStab)

338.18

353.52

3

ParAIP(1)

12.65

parGeo

parBiCGStab

249.75

262.40

7

ParAIP(1)

5.91

parGeo

parBiCGStab

120.34

126.25

15

ParAIP(1)

3.40

parGeo

parBiCGStab

65.41

68.81

Task 9

ExxonMobil’s Geo problem. (geo_32_64_16)

Matrix size: N=244081, NNZ=6218735

- Central difference schemes for Euler gas dynamics (Kurganov-Tadmor difference schemes)
- Application of KT difference scheme to gas pipe simulation problem
- KT scheme for 1D Lagrangian gas dynamics
- KT scheme for multi-dimensional gas dynamics irregular grid simulations

Euler equations

The non-steady-state flow of heat-conductive viscous gas in a constant

section pipe is described by system of differential equations in partial

derivatives [1]:

(1)

Equation of state:

(J/kg К) – universal gas constant,

- molecular weight,

- Poisson index.

- time and space coordinates along pipeline axis.

- density, mass velocity, pressure, specific internal energy and gas temperature

- temperature of environment (ground, water, air)

- pipe diameter,

- pipe radius

- pipeline length (simulation region).

- given function, describing relief along a pipeline.

- resistance coefficient.

- Reynold's number.

- dynamic gas viscousity

- effective rough inner pipe surface

- heat transfer coefficient

Difference scheme

The uniform grid is made along a pipeline section. Density, speed, specific internal energy and temperature of gas are defined in middle of intervals. In addition one accounting interval is attached to external borders of the pipeline, so-called fictitious cells. The known spatial distribution of appropriate values is assign in middle of intervals in the initial time moment.

Second order semi-discrete central difference scheme [2] is used for solution of equations (1):

, where

Difference boundary conditions

It was supposed, that the condition of inflow gas is put on a left border, and the pressure is set on a right border. Linear dependencies of basic values chosen as initial approximation.

1. Inflow condition.

Compositions or are set on inflow. Let's set for definiteness, then temperature we shall express from the equation of state:

Values with index "0" refer to border (node of accounting cell), "R" (right) - to center of the first accounting cell, "L" (left) - to center of the fictitious cell. Let's find values at the center of the fictitious cell. Linear extrapolation of density, impulse and pressure is used for this purpose.

Accordingly, flows on boundary edge will be written down as:

The value of mass velocity is necessary for finding flows.

Small oscillation arise at such assignment order %, as it is visible from the dependence of stationary flow via coordinate, given in figure.

1 variant: linear extrapolation of mass velocity

2 variant: linear extrapolation of impulse

In this case flow is monotone at the left border.

Only pressure is known on outflow.

2. Outflow condition.

Values with index "0" refer to border (node of accounting cell), "L" (left) - to center of the last accounting cell, "R" (right) - to center of the fictitious cell. Let's carry out linear extrapolation of density and impulse. If we set boundary conditions and flows to

then the large entropy trace arises on the right border (up to 1-2%), which scale does not depend on number of grid intervals.

Received equation multiply on , combine with the first equation (1), multiplied on

and subtract the third equation. In result we have,

Subtract last equation from the first equation (1) multiplied on :

Multiply this equation by and combine with the first equation (1) multiplied by :

1. Therefore we have tried to find additional difference boundary conditions being a consequence of the basic boundary condition and equations (1). For this purpose the first equation of system (1) multiply on “u” and subtract from second equation.

, where c is a sound speed

Last equation write down in difference form and take into account, that at , follow , then

The equation below is received using the equation of state of ideal gas

(2)

3. Entropy equation

The law of entropy conservation is used as third difference boundary condition:

The law of entropy conservation in finite differences has the following form:

(3)

Using (2) and (3), one can obtain a system of equations with respect toand

Using the found values and , it is possible to find the other unknowns

Error in temperature at the right boundary decreases with the increasing of grid points number.

Numerical results

- number of difference intervals. Calculations were carried out in interval

Heat exchange with an environment is completely absent in this task: . Calculation results in comparison by exact solutions are given in Fig. 1. Good consent is observed with exact solutions.

The test task given in the report [1], we solve for evaluation of the numerical technique. Gas compressibility, heat exchange with walls, friction on internal surface of a pipe, viscosity essentially influence gas dynamics in these tasks.

Stationary gas flow in heat-insulated pipeline

The stationary task 1 was solved by a method of establishment flow under given boundary conditions.

Pressure

Temperature

Velocity

Mass flow

Fig. 1.Stationary gas flow distribution of pressure, temperature, velocity and

mass flow in in heat-insulated pipeline. Solid line – numerical results,

stars – exact solutions.

In this task the convergence of numerical solution to analytical was analyzed.

Amount of spatial cells is plotted on the blue diagram on X-axis in logarithmic scale, logarithm of norm of a difference numerical solutions and analytical – on Y-axis.

The green colour plots the cubic law: .

It is possible to make a conclusion that the solution converges to analytical with the third order depending on amount of spatial cells.

Literature

1. Anuchin M.G., Dremov V.V. Model for calculate flow of natural gas on gas pipeline. RFNC-VNIITF report, № ПС.99.7432, Snezhinsk, 1999.

2. Kurganov A., Tadmor E. High-resolution central schemes. UCLA Report №16, 1999.

- Problem Statement
- Material Model (by Bychenkov)
- Lagrange Irregular Voronoi Grids
- Central Difference Scheme for Irregular Grids in 2D and 3D Cases

The additional modes of the model

- Porosity accounting
- Material crash accounting
- Welding of the material accounting

- The multidimensional Voronoi Grid is used
- To obtain scheme for multi-dimensional case one should reformulate original 1D KT scheme to make generalization more obvious (see the next 2 slides)

- In the next five slides the results of a test problem numerical modeling are presented
- The initial conditions are: the first region – density 4, specific internal energy – 0.5, velocity – 1; the second region – density – 1, specific internal energy – 1e-6, velocity – 0.
- Boundary conditions: the left boundary – piston with velocity 1, the right one – rigid wall.
- The plots present density and velocity for 10 successive times.

density

velocity

density

velocity

density

velocity

density

velocity

density

velocity

- The numerical results show, that the scheme has inherited its good quality from the Euler parent KT scheme
- The smoothness is minimal (even after more, than 10 shock reflections the number of point at the shock is not more then five)
- Some entropy traces are presented (at the boundaries and at the place, where initially the shock has been located)

- is a set of numbers of cells, which are neighbor cells to the k-th one
- DI, DE are the results of minmod procedure applied to the gradients (in common case six ones), calculated in the Delauneux triangles, including k-th point
- Dr is a result of minmod operation applied to the corresponding laplacians, calculated in the neighbor points

Lagrangian KT Scheme for Irregular Grid