Nonclassified activity of cml rfnc vniitf
Download
1 / 113

Nonclassified activity of CML RFNC-VNIITF - PowerPoint PPT Presentation


  • 160 Views
  • Updated On :

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

loader
I am the owner, or an agent authorized to act on behalf of the owner, of the copyrighted work described.
capcha
Download Presentation

PowerPoint Slideshow about 'Nonclassified activity of CML RFNC-VNIITF' - erno


An Image/Link below is provided (as is) to download presentation

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 - - - - - - - - - - - - - - - - - - - - - - - - - -
Presentation Transcript
Nonclassified activity of cml rfnc vniitf

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 Technical Physics

  • 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.



The main directions of the work
The main directions of the work Technical Physics

  • 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


2 5d mhd code
2.5D MHD Code Technical Physics

  • 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

2½D MHD MAG Code Technical Physics

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 mag code model
The MAG Code Model: Technical Physics

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

r


Details of numeric
Details of Numeric: Technical Physics

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.


1 gas dynamics conditions

Boundary conditions: Technical Physics

1. Gas dynamics 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);


2 conditions for heat conductivity

Boundary conditions: Technical Physics

2. Conditions for heat conductivity:

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.


3 conditions for magnetic field

Boundary conditions: Technical Physics

3. Conditions for magnetic field

• 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 Bis 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.


Algorithms of mesh reconstruction
Algorithms of mesh reconstruction: Technical Physics

  • 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


System of equations in arbitrary moving coordinates i
System of equations in arbitrary moving coordinates I: Technical Physics

Equation of continuity:

Euler equations for velocity components:

x component:

y component:

z component:


System of equations in arbitrary moving coordinates ii
System of equations in arbitrary moving coordinates II: Technical Physics

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=uk, 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.


System of equations in arbitrary moving coordinates iii

r Technical Physics= const. - equation of continuity,

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

u = const. - Euler equation for velocity.

System of equations in arbitrary moving coordinates III:

Equation for A:

Equation for Bz:

Equation for energy:



Diffusive equations system
Diffusive equations system: Technical Physics

Difference sheme for diffusive equations system:


Z pinch simulation
Z-pinch simulation Technical Physics

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.


Hollow gas puff simulation
Hollow gas puff simulation: Technical Physics

Magnetic field,

-component

[10 MGs]

Density [g/cc]


Pinhole image reconstruction throw steady equation of radiation transfer resolving
Pinhole image reconstruction throw steady Technical Physicsequation of radiation transfer resolving

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


Pinhole image reconstruction on simulation of gas puff implosion i
Pinhole Image Reconstruction on Simulation of Technical PhysicsGas Puff Implosion (I)


Pinhole image reconstruction on simulation of gas puff implosion ii
Pinhole Image Reconstruction on Simulation of Technical PhysicsGas Puff Implosion (II)


Laser produced plasma jet expansion into vacuum
Laser produced plasma jet expansion into vacuum Technical Physics

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.


Laser produced plasma jet expansion into background plasma
Laser produced plasma jet expansion into background plasma Technical Physics

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


Experimental setup
Experimental Setup Technical Physics



Plasma channel formation temperature evolution
Plasma Channel Formation. Temperature Evolution. Technical Physics

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.


2d irregular grid for mathematical modeling
2D Irregular Grid for Mathematical Modeling Technical Physics

  • Gas Dynamics

  • Poisson Equation Solver

  • Maxwell Equation Solver

  • Heat Transfer


IGM Code Technical Physics

Oleg V. Diyankov

Sergei S. Kotegov

Vladislav Yu. Pravilnikov

Yuri Yu. Kuznetsov

Aleksey A. Nadolskiy

RFNC – ARITPh

supported by LLNL grant B329117


Overview
Overview Technical Physics

  • 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.


Features benefits
Features & Benefits Technical Physics

  • 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.


Applications
Applications Technical Physics

  • Plasma physics (instability study, laser produced plasma, and so on).

  • High velocity impact.

  • Heat transfer.

  • Electrostatic fields.


Voronoi diagram
Voronoi diagram Technical Physics

Benefits:

  • Local orthogonality

  • Local uniformity




Difference scheme
Difference scheme Technical Physics


Corresponding fluxes
corresponding fluxes Technical Physics


Coordinates
coordinates Technical Physics


Equation of heat conductivity
Equation of heat conductivity Technical Physics

difference scheme:

boundary condition:


Poisson equation
Poisson equation Technical Physics

Integrate over Voronoi cell

where


Electric field strength
Electric field strength Technical Physics

where






Poisson problems degrees)f(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


D ischarger distribution of potential 10 5 v
D degrees)ischarger(distribution of potential, 105 V)


D ischarger distribution of electric field strength 10 5 v cm
D degrees)ischarger(distribution of electric field strength, 105 V/cm)


Specifications
Specifications degrees)

  • 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


  • · Linear Systems Solving Introduction

  • ·Iterative Methods

  • ·Preconditioners

  • ·The PLS Code


Introduction
Introduction Linear Systems Solving

  • 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 equations

, 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 equations

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 equations

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 equations

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 equations

(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 equations

  • 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 equations

  • 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 equations

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 equations

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: equations

  • 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.


Numerical results
Numerical results equations

  • 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.


N equationsproc

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.


N equationsproc

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


N equationsproc

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


N equationsproc

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


N equationsproc

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


N equationsproc

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


N equationsproc

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


N equationsproc

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


N equationsproc

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


Difference schemes for hyperbolic systems treatment
Difference Schemes for Hyperbolic Systems Treatment equations

  • 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 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. equations

- 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 equations

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 equations

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: equations

The value of mass velocity is necessary for finding flows.


Small oscillation arise at such assignment equationsorder %, 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. equations

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 equations, 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 account, that at

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 account, that at 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 account, that at

- 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 account, that at

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 analytical was analyzed.

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.


3d elastic plastic modeling of processes of ceramics formation
3D Elastic-Plastic Modeling of Processes of Ceramics Formation

  • 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 Formation

  • Porosity accounting

  • Material crash accounting

  • Welding of the material accounting







Difference scheme for multidimensional irregular grid
Difference scheme for multidimensional irregular grid Formation

  • 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)


Difference scheme for 1D Lagrangian Gas Dynamics (Variation for generalization for multidimensional irregular grid)


Difference scheme for 1D Lagrangian Gas Dynamics (Variation for generalization for multidimensional irregular grid)


A test problem for the lagrangian difference scheme
A test problem for the lagrangian difference scheme for generalization for multidimensional irregular grid)

  • 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 for generalization for multidimensional irregular grid)

velocity


density for generalization for multidimensional irregular grid)

velocity


density for generalization for multidimensional irregular grid)

velocity


density for generalization for multidimensional irregular grid)

velocity


density for generalization for multidimensional irregular grid)

velocity


Numerical results discussion
Numerical results discussion for generalization for multidimensional irregular grid)

  • 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)


Explanation of some details of the multidimensional scheme presented on the next slide
Explanation of some details of the multidimensional scheme, presented on the next slide.

  • 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 presented on the next slide.


ad