MA375 - Rice U - Computational and Applied Mathematics

Download Report

Transcript MA375 - Rice U - Computational and Applied Mathematics

MA557/MA578/CS557
Lecture 22
Spring 2003
Prof. Tim Warburton
[email protected]
1
Interpolation on the Triangle
• Recall we are considering a two-dimensional domain.
• We assume that a triangulation of the domain is given.
• In each triangle we are going to create a polynomial
approximation of the solution to a PDE.
• We discussed briefly an orthogonal basis for the triangle
(which you are currently verifying is indeed actually an
orthogonal set of functions).
• We construct a generalized Vandermonde basis using
this basis and a set of nodes.
2
Reference Triangle
• The following will be our basic triangles:
(-1,1)
s
r
(-1,-1)
(1,-1)
• All straight sided triangles are the image of this
triangle under the map:
1
2
3





v
v
v
 x
 r  s  x  1 r  x
 1 s  x 
 y     2   v1    2   v 2    2   v3 

 y  
 y  
 y 
 
3
Reference Triangle
v , v 
3
x
• The following will be our basic triangles:
(-1,1)
s
r
(-1,-1)
3
y
v , v 
1
x
1
y
v , v 
2
x
(1,-1)
• All straight sided triangles are the image of this
triangle under the map:
1
2
3





v
v
v
 x
 r  s  x  1 r  x
 1 s  x 
 y     2   v1    2   v 2    2   v3 

 y  
 y  
 y 
 
4
2
y
• Given a set of nodes lying in the triangle we use V to
construct an interpolating polynomial for a function who’s
values we know at the nodes:
i  p j  p i
f  r, s   
i 0
ˆ

r
,
s
f


 ij 
 ij 
j 0
• The interpolation condition yields:
i  p j  p i
f  rn , sn   
i 0
    r , s  fˆ 
j 0
ij
n
n
ij
i  p j  p i

i 0

m M
 V   fˆ 
j 0
V
m 1
n ij
ij
fˆ where we have chosen an ordering of the modes
nm m
  p  1 p  2  
M 

2


5
Differentiation
• Suppose we wish to find the derivative of a p’th order
polynomial
f  P p T  for a triangle T
• First we note that the approximation becomes equality:
i  p j  p i
f  r, s  
  r , s  fˆ

i 0
j 0
 ij 
 ij 
• And interpolation allows us to find the PKDO coefficients:
f  rn , sn  
m M
V
fˆ
nm m
m 1
• So differentiation requires us to compute:
i  p j  p i 
f
 ij 
 r, s    
 r , s  fˆij 
r
r
i 0 j 0
i  p j  p i 
f
 ij 
r
,
s

  
 r , s  fˆij 
s
s
i 0 j 0
6
Differentiation cont
• So we need to be able to compute:
 ij 
r
and
 ij 
s
 r, s 
• Recall the definition of the basis functions:
  1  s  2 n1,0
0,0  2(1  r )
 nm  r , s   Pn 
 1 
s
 Pm
 (1  s)
 2 
n
• R-derivative:
 nm
 2  dP  2(1  r )   1  s  2 n1,0
 1 
 r, s   
s

 Pm

r
 1  s  dx  (1  s)
 2 
0,0
n
n
7
Quick Jacobi Polynomial Identity
• We will make extensive use of the following:
dPn ,
 n      1   1, 1
 x  
 x
 Pn1
dx
2


8
r-Derivative
• Ok we need to calculate:
 nm
 2  dP  2(1  r )   1  s  2 n1,0
 1 
 r, s   
s

 Pm

r
 1  s  dx  (1  s)
 2 
n
0,0
n
 nm 
 2  n  1  1,1  2(1  r )   1  s  2 n1,0
 1 
 r, s   
s n  1

 Pn1 
 Pm
r
 1  s  2 
 (1  s)
 2 
n
 0 when n  0 (since P0 ,  x   1)
• We can compute these using the definition of the Jacobi
polynomials.
• Watch out for s=1 (top vertex) – the r-derivative of all the
basis is functions is zero at r=1,s=-1
9
s-Derivative
We use the chain and product rule to obtain:
 2 1  r  dP 0,0  2(1  r )    1  s n 2 n1,0
n
 1  
 r, s   
s
 Pm

2
s
   2 
 1  s  dx  (1  s )
 nm 
  n  1  s 
0,0  2(1  r )
 Pn 
 1  

(1

s
)
2
2


  
n 1
 2 n1,0
s
 Pm

2 n 1,0




dP
2(1

r
)
1

s


0 ,0
m
 Pn 
 1 
 s 
 
 (1  s )
  2   dx

n
10
s-Derivative
From which:
 2 1  r   n  1  1,1  2(1  r )    1  s  n 2 n1,0
 1  
 r, s   
s
 Pn1 
 Pm
2 
s
 (1  s )
   2 
 1  s   2 
 nm 
  n  1  s 
0,0  2(1  r )
 Pn 
 1   

(1

s
)
2
2


  
n 1
 2 n1,0
s
 Pm

  1  s   m  2n  2  2 n2,1 
0,0  2(1  r )
 Pn 
 1 
 
 Pm1  s  
2


 (1  s )
  2  
n
11
Special Cases
 0 m 
 1,0
 r , s   Pm  s 
s
s
Don’t worry about all those denominators having (1-s)
since the functions are just polynomials and not
singular functions…
12
Recap
 nm 
 2  n  1  1,1  2(1  r )   1  s  2 n1,0
 1 
 r, s   
s n  1

 Pn1 
 Pm
r
 1  s  2 
 (1  s)
 2 
n
 0 when n  0 (since P0 ,  x   1)
 2  n  1  1,1  2(1  r )    1  s  n 2 n1,0
 1  
 r, s   
s
 Pn1 
 Pm
2 
s
 (1  s )
   2 
 1  s   2 
 nm 
  n  1  s 
0,0  2(1  r )
 Pn 
 1  

(1

s
)
2
2



 
n 1
 2 n1,0
s
 Pm

 2(1  r )   1  s   m  2n  2  2 n2,1 
 Pn0,0 
 1 
 
 Pm1  s  
2


 (1  s )
  2  
n
13
Derivative matrices
• Given data at M=(p+1)(p+2)/2 points we can directly r
and s derivatives with:
fˆm 
m M
1
V
 
m1
mj
f  rj , s j 
m M
m
f
r ˆ
r
ˆ
ˆ
 rn , sn    Dnm f m where Dnm 
 rn , sn 
r
r
m1
m M
m
f
s ˆ
s
ˆ
ˆ
r
,
s

D
f
where
D

 n n   nm m
 rn , sn 
nm
s
s
m1
14
One-Stage Differentiation
• Given a vector of values of f at a set of nodes we can
obtain a vector of the r and s derivatives at the nodes
by:
f
ˆ r V 1
 Dr f where Dr  D
r
f
s
s
s 1
ˆ
 D f where D  D V
s
15
Matlab Scripts
• After class on Friday I will post code which computes
the derivative matrices Dr and Ds for an arbitrary set
of M nodes.
16
Inner Product Matrices
• Recall we need to compute:
 f , g T  
fgdV
T
  x, y 
   f  r, s  g  r, s 
drds
  r, s 
1 1
1 s
• It is not obvious how to do this given the value of f
and g at a set of M points.
• Same old trick – construct the PKDO coefficients and
use the orthogonality relationship..
17
Inner-Product
 f , g T
  x, y  1  s

f  r , s  g  r , s  drds


  r , s  1 1
Since the r,s->x,y is linear
  x, y  1  s  m  M ˆ
 nM

ˆ

f

r
,
s
g

r
,
s
  m m     n n    drds


  r , s  1 1  m1
 n1

  x, y  1  s  n  M m  M ˆ

ˆ

   f mm  r , s g nn  r , s   drds


  r , s  1 1  n1 m1

1 s
  x, y  n  M m  M ˆ

f m gˆ n   m  r , s n  r , s  drds


  r , s  n1 m1
1 1
  x, y  n  M m  M ˆ

2
 2 
ˆ

f
g

C
where
C

  m n nm n


 ij 
  r , s  n1 m1
2
i

1
2
i

2
j

2



  x, y  n  M ˆ

f n gˆ nCn

  r , s  n1
18
Nodal Mass Matrix
• Suppose we know the value of f and g at M points
then we can compute the PKDO coefficients using
the generalized Vandermonde matrix.
• We can then integrate by the previous operations
(last slide).
• All these operations can be concatenated into one
mass matrix.
19
Mass-Matrix
 f , g T
  x, y  n  M ˆ

f n gˆ nCn

  r , s  n1



  x, y 
  r, s 
  x, y 
  r, s 
  x, y 
  r, s 
1
V
f
C
V
   g
t
1
f V
t
t

1 t
f Mg
where C is a diagonal matrix
C  V 1  g
where M   V
 CV 
1 t
1
20
Nodal Mass Matrix
• Setting
• Where:
f  hn  r , s 
g  hm  r , s 
hn  rm , sm    nm 1  n, m  M ,
 hn , hm T
p  1 p  2 

where M 
2
 htn Mh m

• Then:

  x, y  k  M
j M
h r , s M


  r, s 
n
k 1
j 1
  x, y  k  M
j M



  r, s 
k 1

  x, y 
  r, s 
nj
j
j
h  rk , sk 
jk m
M jk  km
j 1
M nm
21
DG Matrices
 g 
 f, 
 x T
• Recall we need to compute:
• We use the coordinate change, chain rule, linearity of the
map T->That (reference triangle) and finally the identities we
just found:
 g   r g   s g 
 f,   f,
  f ,


x

x

r

x

s

T 
T 
T
   x, y  r g     x, y  s g 
 f,
  f ,

   r , s  x r Tˆ    r , s  x s Tˆ
  x, y  r  g    x, y  s  g 

 f,  
 f, 
  r , s  x  r Tˆ   r , s  x  s Tˆ

  x, y  r
  r , s  x
f MD g 
t
r
  x, y  s
  r , s  x
f t MD s g
22
Summary
Set:   r , s   hn  r , s  where hn is the n’th Lagrange interpolant
defined as the p’th order polynomial in r,s for which:
hn  rm , sm    nm 1  n, m  M ,
p  1 p  2 

where M 
2
   x, y  r 
   x, y  s 
 hm 
r
s
h
,

MD

MD








 n

nm
nm
x T    r , s  x 

   r , s  x 
   x, y  r 
   x, y  s 
 hm 
r
s
h
,

MD

MD








 n y 
nm
nm
   r , s  y 

T    r , s  y 
23
Progress With The DG Scheme For Advection
• The DG scheme now looks like:
Find Cm  t  , m  1,..., M  such that:
m M
 h ,h 
n
m 1
m T

 mM   hm 
dCm mM   hm 
  u  hn ,
Cm .....
 Cm    v  hn ,

dt
x T  m1  
y T 
m 1  
     
  hn , 
  2
 


 C    0



T
where    u  n  ,
C   C   C 
• Where we set:   r , s   hn  r , s  and we can now compute
everything in the first three inner-product.
• Next time we will discuss how to compute everything in the
surface inner-product.
24
Summary of Matrices
   x, y  r 
   x, y  s 
 hm 
r
s
h
,

MD

MD






 n

nm
nm
x T    r , s  x 

   r , s  x 
   x, y  r 
   x, y  s 
 hm 
r
s
h
,

MD

MD








 n y 
nm
nm

T    r , s  y 
   r , s  y 
 hn , hm T
   x, y  

M
   r , s  
r r s s   x, y 
, , , ,
Recall the factors: x y x y   r , s 
over a straightsided triangle.
are constant
25
Next Lecture
• We will use a better set of nodes (improve the
condition number of the Vandermonde matrix).
• Time permitting we will prove consistency.
26