Finite Volume Time Domain Technique (FVTD) Computations

Download Report

Transcript Finite Volume Time Domain Technique (FVTD) Computations

Finite Volume Time Domain (FVTD)
Computations for Electromagnetic
Scattering
Prof. A. Chatterjee
Department Of Aerospace Engineering
IIT Bombay.
Outline Of Presentation
Part I : Introduction to the finite volume time domain technique
Maxwell Equation in Conservative form
Finite volume method for conservation laws
Numerical Schemes
Volume grid generation
Boundary Conditions
Validation
PEC sphere
Almond
Ogive
lossy sphere
lossy cone sphere
Part II: Applications
RCS prediction of low observable aircraft Configurations
Intake
B2
F117
Maxwell Equation (differential form)
Maxwell’s curl Equations with losses:
B
   E  H
t
D
   H  J  E
i
t
Constitutive Relations
D   E
B   H
   
   
 r     j 
 r     j 
Maxwell Equation in Conservative form
3D Maxwell’s Equations in Conservative form:
ut  f (u) x  g(u) y  h(u) z  s
where,
 Bx 
  Dy /      H x 
 0 
Dz /   


 

 




 By 
 Dx /      H y 
  Dz /   
 0 
B 
 0    H 
 D / 




D
/

y
z
z 


; s  

; g   x ; h  
u
;f 
 Dx 
 B y /      J ix  E x 
 0 
  Bz /   


 






0
 J iy  E y 


B
/

x

 Dy 


 Bz /   




  B /  
 0    J  E 
D 

B
/



y
x



  iz
 z
z
Numerical Technique
Finite Volume Time Domain (FVTD) technique.
Higher Order Characteristic based technique for spatial discretization
based on the Essentially Non-Oscillatory (ENO) method.
Multi-stage Runge-Kutta time integration.
Numerical Formulation
Maxwell’s Equations (Operator form)
L(u )  s
Decomposition of Total Field
L(u  u )  s  s
i
s
i
L(u )  s  S
s
s
i
S   L(u )  s
i
i
i
s
Finite Volume Framework
Conservative form can be written as
u  f (u ) x  g(u ) y  h(u ) z  s  S
s
t
s
s
s
s
i
where
S i  u it  f (u i ) x  g(u i ) y  h(u i ) z  s i
integrating over an arbitrary control volume,
  u s dV
v
t
   [F (u s )]dV
v
  ( s s  s i )dV 
v
  u i dV
v
t
   [F (u i )]dV
v
Finite Volume Framework ….contd
application of divergence theorem gives,
  u s dV
v
t
  [F(u s )  F(u i )]  nˆdS   ( s s  s i  uti )dV
s
v
3D domain divided into hexahedral cells
Cell centred formulation
discretized form for jth cell,
Vj
~s
du
j
dt

M
 ({[F (u s )  F(u i )]  nˆ S }m ) j  V j (~s js  ~s ji 
m 1
Higher order characteristic based technique
Runge Kutta time stepping
~i
du
j
dt
)
Boundary Conditions
On Perfectly Conducting (PEC) surface
Total tangential electric field,
n E  0.
Total normal magnetic field,
n  B  0.
Far-field boundaries
Characteristic boundary conditions (zero scattered field)
Numerical Boundary Treatment
Surface of the body is a Perfect Electric Conductor (PEC)
Tangential component of electric field zero at surface, i.e, n E  0.
In scattered field formulation, implemented as n  Escat  n  Einc (Einc known)
Similarly the boundary condition for the magnetic field, n  B  0.
Field values in the ghost cells computed by extrapolating the scattered field values from the interior
Normal components of the electric fields and tangential components of the magnetic fields in the
ghost cells taken identical to those on the surface of the conductor
Boundary Conditions ….cont.
Numerical Boundary Treatment
In the far field, characteristic based boundary treatments are applied
Scattered field is taken as zero
Fluxes are decomposed along the characteristic directions normal to the
cell faces; for the ghost cells outside the domain, incoming characteristic
fluxes for the scattered field are taken to be zero
Methodology
Time domain computations for sinusoidal steady state
Complex field in frequency domain from time history of solution using
Fourier Transform
Finite Volume Formulation for
Conservation Laws
Introduction
Aim: To solve a set of governing equations describing a set of
conservation laws in the integral form in a specified domain with
prescribed boundary conditions
Equation of Conservation Laws
e.g. In Fluid Mechanics
Mass conservation
Momentum conservation
Energy conservation
+ Equation of state as a closure equation
e.g. In Electromagnetics
Maxwell’s Curl Equations
e.g. In Magnetohydrodynamics
Navier Stokes Equations & Maxwell’s Equations
Conservation Equation
Consider generic conservation equation for a conserved variable u, and
assume that the corresponding flux vector known
ut  f (u) x  g(u) y  h(u) z  s
nˆ
Integrating over an arbitrary volume
  udV
v
t
  [F(u)]  nˆ dS   ( s )dV
s
v
where,
F  f(u)iˆ  g(u) ˆj  h(u)kˆ
Conservation Equation for 
The conservation equation for  in a general form:
 

dV   V  S dS   q dV

t V
S
V
Rate of change
of  in a control
volume
Convective
flux across
the surface
Here :  = 1 Mass conservation equation
= u,v,w Momentum conservation equation
= Energy equation
Source or
sink term
Finite Volume Method for Conservation Laws
Methodology
Starting point: Integral form of conservation equation.
Domain covered by finite number of contiguous control volumes (CV)
Conservation equation is applied to each CV
Computational node : Center of each CV, where the Variable values are
calculated.
Interpolation is used to calculate variable values of CV surfaces in terms of
nodal (CV center) values.
Results in an algebraic equation for each CV per equation
Computational Domain for Finite Volume Method
Computational domain
divided into finite
control Volumes
Surface Grid
Control Volume
Volume Grid
Finite Volume Discretizations
Geometry Creation
Domain Discretization: Domain is subdivided into a finite number of
small control volumes (CV) by a grid which defines the control volume
boundaries
Grid Terminology
nˆ
Typical CV with the notations
Node-based finite volume scheme: u stored at vertex
Cell-based finite volume scheme: u stored at cell centroid
Finite Volume Method … contd.
The generic conservation equation for a conserved variable u, reproduced
here is applied to each CV
  udV
v
t
  [F(u)]  nˆ dS   ( s )dV
s
v
Net flux through CV boundary faces =  f  dS   f num  S
k
S
 

  V  


f = component of the convective
direction normal to the CV face
flux vector in the
Note: CVs should not overlap
Domain volume should be equal to the sum of volumes of all the CVs
Each CV face is unique to two CVs which lie on either side of it
Finite Volume Method … contd.
To evaluate surface integral  f  dS   f num  S exactly one
S
k
should know integrand f everywhere on “S”
The above information is not available, only nodal (CV center) values
of ‘u’ are calculated. Approximation must be introduced
Integral is approximated in terms of the variable values at one or more
locations on the cell face.
Cell face values are approximated in terms of nodal CV center values
Finite Volume Method … contd.
Global Conservation:
Integral conservation equation is applied to each CV
Sum equations for all CVs
Global conservation is obtained since surface integral over CV faces
cancel out
Advantage
Simplest method
Dis-advantage
Three levels of approximations
Complex geometries
• Interpolation
Scheme is conservative
• Differentiation
Most popular amongst CFD workers
• Integration
( > 95%)
Numerical Schemes
Differ by how numerical flux fnum is evaluated at cell face and by how
time integration is performed
Differ in spatial and temporal order of accuracy
Upwind (Characteristic based) Schemes:
Flux Splitting Schemes, Riemann/Godunov Solvers
Steger Warming, Van Leer flux vector splitting
Roe type solvers
Central Difference based Schemes:
Lax-Wendroff Scheme
Jameson’s scheme
Time Stepping
Space and Time discretization combined
e.g. Lax Wendroff Scheme
Taylor Series expansion in Time
Space and Time seperated
Set of ODE’s obtained after space discretization
du
 R (u )  0
dt
R (u ) 
1
 f  S
V
March in time with Runge-Kutta method
ENO (Essentially Non-Oscillatory)Scheme
Higher order accuracy
Non-oscillatory resolution of discontinuities (adaptive stenciling)
Models discontinuous changes in material properties
Numerical flux at right-hand face of jth cell with rth-order accurate ENO
scheme
r 1
f j 1 / 2  
l 0
r
k ,l
f j r 1 k l
The smoothest possible stencil
Sk  ( x j k r 1 , x j k r 2 ,.....,x j k )
Validation
for Standard Test Case
Validation
Metallic Sphere:
Volume Grid – O-O Topology, Single Block (50×45×20 cells)
Frequency for analysis = 0.09 GHz ( Electric Size = 1.4660 )
Volume Grid Discretization (O-O Topology)
Validation ….contd
Metallic Sphere:
E-plane
h-plane
Bistatic RCS (dB) for Metallic Sphere at 0.09 GHz
Validation - Metallic Ogive
Metallic Ogive:
Electro-Magnetic Code Consortium (EMCC) benchmark target
Geometry – Half angle 22.62 degrees, length 5”, thickness 1”
Volume Grid - O-H topology, single block (40×40×138 cells)
Monostatic RCS at 1.18 GHz (Electric size 2π) for vertical polarization
Excellent agreement with experimental results
Maximum RCS of approx -20 dBsm = 10-2 m2
Minimum RCS of approx -60 dBsm = 10-6 m2
Validation ….contd
Metallic Ogive:
Surface Grid
Rendered Ogive
Validation ….contd
Metallic Ogive:
Ogive Volume Grid Cross-Section
Ogive – Volume Grid and Surface Currents
Validation ….contd
Metallic Ogive:
Ogive Monostatic RCS Plot (1.18 GHz, VV Polarization)
Validation ….contd
Almond:
Electro-Magnetic Code Consortium (EMCC) benchmark target
Geometry – length 9.936”
Volume Grid - O-O topology, single block (15×121×125 cells)
Monostatic RCS at 1.19 GHz (Electric size 2π) for vertical polarization
Excellent agreement with experimental results
Validation ….contd
Almond:
Surface Grid
Rendered Almond
Validation ….contd
Almond:
Angle of incidence = 0 deg.
Angle of incidence = 90 deg.
Almond – Surface Currents
Validation ….contd
Almond:
Almond Monostatic RCS Plot (1.19 GHz, VV Polarization)
Validation ….contd
PEC Sphere with Non-Lossy coating:
PEC ka = 2.6858
Coating (t / λ) = 0.05, ka1 = 3.0, ε‘ = 3.0 and 4.0, μ‘ = 1.0
Volume Grid O-O Topology, Single Block (64×45×32 cells)
Bistatic RCS for Sphere with Discontinuous Nonlossy Coating
(Backscatter at 180 degree)
Validation ….contd
PEC Sphere with Lossy dielectric coating:
PEC ka = 1.5
Coating (t / λ) = 0.05, εr = 3.0 – j4.0, μr = 5.0 – j6.0
Volume Grid O-O Topology, Single Block (64×48×32 cells)
For different orders of accuracy
For different discretization
Monostatic RCS Sphere with Lossy Coating
Validation ….contd
PEC Cone Sphere with Lossy Coating:
Geometry : vertex angle 90 degree, sphere diameter 0.955 λ
Coating (t / λ) = 0.01, εr = 3.0 – j4.0, μr = 5.0 – j6.0
Volume Grid O-O Topology, Single Block (80×40×38 cells)
E-plane
h-plane
Bistatic RCS for cone-sphere with lossy Coating
(Backscatter at 180 degree)
Radar Cross Section of Low
Observable Aircraft Configurations
Applications: Industrial Problems
Industrial Problems
RCS Analysis of Engine intake configurations
B2 “Advanced Technology Bomber”
F-117 “Nighthawk”
RCS of Low Observable Aircraft Configurations
Introduction
RCS, low-observability, air intake configurations
Intake geometries and grid generation
Results
RCS Analysis of Engine Intake configurations
Introduction:
Low-observability: low back-scatter for near-axial incident illumination as well as low
returns over a broad angular region
Low-observables characterized by greater contribution of traveling and creeping waves to the
Radar Cross-Section (RCS)
Definition of RCS: the area of an isotropic reflector returning the same power per solid
angle as the given body. At far field, it is proportional to the ratio of the power received from
the target to the power incident on the target.
Engine Intake Configurations studied
B2 “Advanced Technology Bomber” and F-117 “Nighthawk” chosen for present study:
Low-observable (stealthy) aircraft configurations (RCS of -40 and -25 dBsm respectively)
Fine geometric details not available due to military sensitivity
RCS of some military aircrafts
Aircraft
RCS
RCS
RCS
[dBsm]
[m2]
[ft2]
F-15 Eagle
26
405
4,358
F-4 Phantom II
20
100
1,076
B-52 Stratofortress
20
99.5
1,071
Su-27
12
15
161.4
B-1A
10
10
107.6
F-16 Fighting Falcon
7
5
53.82
B-1B Lancer
0.09
1.02
10.98
F-18E/F Super Hornet
0
1
10.76
BGM-109 Tomahawk
-13
0.05
0.538
SR-71 Blackbird
-18.5
0.014
0.15
F-22 Raptor
-22
0.0065
0.07
F-117 Nighthawk
-25
0.003
0.03
B-2 Spirit
-40
0.0001
0.01
Boeing Bird of Prey
-70
0.0000001
0.000008
Estimation of RCS
Computation of first time-domain numerical solution proposed by Yee in 1966 (FDTD,
second order in space & time); followed by other FDTD-based algorithms
FVTD-based schemes adopted to handle more complicated geometries
Hyperbolicity of the Maxwell's equations in their conservative form exploited by
characteristic-based algorithms
Drawbacks:


Requirement of large computational resources (both processor speeds and memory)
Lack of theoretical estimates on grid-fineness and minimum distance to the far-field
Similar studies:

RCS of VFY218 (Conceptual aircraft): approx. 15 dBsm @ 100 MHz, nose-on
incidence (monostatic) Not a low-observable

RCS of F117 (Stealth Fighter): approx. -20 dBsm @ 215.38 MHz, nose-on incidence
(monostatic)
RCS Calculation
Above algorithm used to compute the total fields at the surface of the PEC.
Equivalent surface currents are given by
ˆ H
Js  n
ˆ E
Ms   n
(equivalent electric current)
(equivalent magnetic current)
Far-zone transform used to calculate the scattered fields Esc and Hsc at
infinite distance and RCS is calculated as
| E sc |2
  4R
| Einc |2
2
where R is taken to be a sufficiently large number (100,000)
RCS ANALYSIS OF ENGINE INTAKE
CONFIGURATIONS
RCS And Air Intake configurations
Width of air intake duct vis-à-vis wavelength of electromagnetic radiation
Sr-71 air intake
Grilled air intake of f117
Position / profile of intake duct
s-shaped wing mounted intake duct of b-2
Wing blended intake duct of yf-23
SR-71
Intake centre
Body
CENTRE
BODY
ANNULAR
REGION
DUCT WALL
F-117 AIRCRAFT
GRILL
ENG.
FACE
SOFT
HARD
FACE
FACE
INTAKE GRILL
B-2 BOMBER
S-shaped duct
with coated wall
YF-23 ADVANCED
TACTICAL FIGHTER
RCS Analysis of Engine Intake configurations
Representative models of engine intake cavities considered are
Straight cylindrical cavity, open from one end and the other end
Terminated by a flat plate
Terminated by a hub and a flat plate
Terminated by a hub with a set of straight rotor blades and a flat plate
 Hub of cylindrical shape
 Hub of sphere-cylindrical shape
Multiblock representation of a straight cylindrical
cavity terminated by a plate
2λ
20 λ
j
k
4λ
i
OPEN END
ƒ = 15 GHz
λ = 20 mm
PLATED END
es = 4 π
Contd ..
Bistatic RCS: Straight Cylindrical Cavity with a
Flat Plate Termination
Monostatic RCS: Straight Cylindrical Cavity with a
Flat Plate Termination
Front view
Φ
Rear view (angle of irradiation φ=28°)
Straight cylindrical cavity with a hub &
plate termination
2λ
k
i
4λ
λ
λ
PLATED END
HUB (SPHERE-CYLINDER)
j
OPEN END
λ = 20 mm
es = 4 π
Bistatic RCS: Straight Cylindrical Cavity with different
terminations and that of a Hub
FRONT VIEW
Straight cylindrical cavity with a hub blade &
plate termination
6λ
45°
k
11
i
i
3λ
6λ
10
5
9
6
4λ
8
7
2.5°
j
OPEN END
4
BLADE
ƒ=6 GHz
HUB
λ=50 mm
PLATED END
es = 6π
VIEW-OPEN
END
One
block
VOL. DISCRETIZATION OF BLOCKS 4-11
(REPRESENTING REGION BETWEEN
THE BLADES)
Bistatic RCS: Straight Cylindrical Cavity with a Hub, Blades
and Plate termination
8GHz, Vertical-Polarization
Monostatic RCS: Straight Cylindrical Cavity with a Hub, Blades and
Plate termination
8GHz, Horizontal-Polarization
Monostatic RCS: Straight Cylindrical Cavity with a Hub, Blades
and Plate termination
6 GHz, Horizontal Polarization, Φ=0°
Surface Current Distribution on Straight Cylindrical Cavity with Hub,
Blades and Plate terminations
8 GHz, Vertical Polarization, Φ=0°
Surface Current Distribution on Straight Cylindrical Cavity with Hub,
Blades and Plate terminations
Straight cylindrical cavity with a hub blade & plate
termination
6λ
45°
k
11
i
4
10
5
9
6
3λ
6λ
8
4.5 λ
2.5°
j
OPEN END
7
BLADE HUB PLATED END
ƒ=6 GHz
λ=50 mm
es = 6π
VIEW-OPEN END
Surface discretization of the complete domain
2
12
1
13
15
14
3
4 -11
Resolution = 20 grid pts/λ
8GHz Vertical Polarization
Comparative Monostatic RCS: Intake cavities with
Cylindrical and Sphere-cylindrical Hub configurations
Side view of the cavity (φ =0)
Front view of the cavity
B2 – Advanced Technology Bomber
Geometry and Grid Generation:
Geometrical details taken from public literature - 3-view diagram from Jane's All the World's
Aircraft
Lack of fine details due to military sensitivity; photographs from the internet used to reconstruct
the geometry
Basic dimensions: 20.9 m length, 5.1 m height and 52.43 m wingspan
Shape of aircraft contributes to its low-observability: blended wing-body, buried engines, lack
of external tanks or weapons, etc
Surface grid: 25 2-D blocks, 45,840 nodes, formed by piece-wise bilinear surfaces
Airfoil section used for the wing: E180 (used in a radio-controlled model of the aircraft)
Average cell dimensions: 0.9 mm in longitudinal direction and 1.7 mm in span-wise direction
Volume grid: 52 blocks, approx. 1.5 million cells, nodes clustered near the surface
Wingspan of the grid: 293.9 mm
Minimum wavelength allowable on present grid (satisfying Nyquist's sampling criterion) is 1.8
mm (electric size 1025.4), corresponding to a wavelength of 0.321 m (frequency of approx 1
GHz) for the full scale aircraft
B2 – Advanced Technology Bomber
Rendered Image of B2 Surface
Grid
Image of the B2 in flight
B2 – Advanced Technology Bomber
B2 Surface Grid (close view):
B2 – Advanced Technology Bomber
Surface Currents on B2 – Nose-on incidence @ 300 MHz,
VV polarization
B2 – Advanced Technology Bomber
Surface Currents on B2 – Broadside incidence @ 300 MHz,
VV polarization
B2 – Advanced Technology Bomber
Surface Currents on B2 – Nose-on incidence @ 300 MHz,
HH polarization
B2 – Advanced Technology Bomber
Surface Currents on B2 – Broadside incidence @ 300 MHz,
HH polarization
B2 – Advanced Technology Bomber
Bi-static RCS Plot @ 50 MHz (VV)
B2 – Advanced Technology Bomber
Bi-static RCS Plot @ 300 MHz (VV)
B2 – Advanced Technology Bomber
Bi-static RCS Plot @ 300 MHz (HH)
B2 – Advanced Technology Bomber
Monostatic RCS Plot @ 300 MHz
F-117 “Nighthawk”
Rendered Image of F-117 Surface
Grid
F-117 “Nighthawk”
F-117 Surface Currents:
Surface Currents at 68 degrees angle of
incidence (HH polarization)
Surface Currents at 90 degrees angle of
incidence (VV polarization)
F-117 “Nighthawk”
Bi-static RCS Plot @ 280 MHz (VV)
Acknowledgement
Anuj Shrimal
Narendra Deore
Manoj Vaghela
Debojyoti Ghosh
Wing Cdr. A. Bhattacharya
Sandip Jadhav
IITZeus Grid Generator
Prof. G.R.Shevare and his Team
Thank you