Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
import numpy as np
import numpy.linalg as la # noqa: F401
import pyopencl as cl # noqa: F401
import pyopencl.array # noqa
import pyopencl.tools # noqa
import pyopencl.clrandom # noqa
import loopy as lp # noqa
ndim = 3
nvars = 3
gamma = 1.4
def rho(state):
return state[0]
def velocity(state):
return state[1:ndim+1]/rho(state)
def energy(state):
return state[ndim+1]
def kinetic_energy(state):
vel = velocity(state)
return 0.5*rho(state)*np.dot(vel, vel)
def internal_energy(state):
return energy(state) - kinetic_energy(state)
def pressure(state):
return (gamma - 1)*internal_energy(state)
def sound_speed(state):
return np.sqrt(gamma*pressure(state)/rho(state))
def flux(state):
vel = velocity(state)
p = pressure(state)
result = np.outer(state, vel)
for i in range(ndim):
result[i+1,i] += p
result[ndim+1,i] += vel[i]*p
return result