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
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
|
C This subroutine sets up Minkowski spacetime with a gague wave.
C
c $Header$
C
C Author: unknown
C Copyright/License: unknown
C
#include "cctk.h"
#include "cctk_Arguments.h"
#include "cctk_Parameters.h"
#define Pi (4 * atan(1.d0))
subroutine Exact__Minkowski_gauge_wave(
$ x, y, z, t,
$ gdtt, gdtx, gdty, gdtz,
$ gdxx, gdyy, gdzz, gdxy, gdyz, gdzx,
$ gutt, gutx, guty, gutz,
$ guxx, guyy, guzz, guxy, guyz, guzx,
$ psi, Tmunu_flag)
implicit none
DECLARE_CCTK_PARAMETERS
DECLARE_CCTK_FUNCTIONS
c input arguments
CCTK_REAL x, y, z, t
c output arguments
CCTK_REAL gdtt, gdtx, gdty, gdtz,
$ gdxx, gdyy, gdzz, gdxy, gdyz, gdzx,
$ gutt, gutx, guty, gutz,
$ guxx, guyy, guzz, guxy, guyz, guzx
CCTK_REAL psi
LOGICAL Tmunu_flag
c local static variables
logical firstcall
CCTK_REAL a, o, d, fs
data firstcall /.true./
save firstcall, a, o, d, fs
c local variables
CCTK_REAL H
c constants
CCTK_REAL zero,half,one
parameter (zero = 0.0d0, half=0.5d0, one=1.0d0)
C This is a vacuum spacetime with no cosmological constant
Tmunu_flag = .false.
C Get parameters of the exact solution.
if (firstcall) then
a = Minkowski_gauge_wave__amplitude
o = Minkowski_gauge_wave__omega
d = Minkowski_gauge_wave__lambda
fs = Minkowski_gauge_wave__phase
firstcall = .false.
end if
C How should the wave look like.
if (CCTK_Equals(Minkowski_gauge_wave__what_fn,"sin").eq.1) then
d = Minkowski_gauge_wave__lambda * half / Pi
if (Minkowski_gauge_wave__diagonal.ne.0) then
H = one - a * sin((x-y)/d - o*t/d)
else
H = one - a * sin((x-o*t)/d - fs)
end if
elseif (CCTK_Equals(Minkowski_gauge_wave__what_fn,"expsin").eq.1) then
H = exp(a*sin(x/d)*cos(t/d))
elseif (CCTK_Equals(Minkowski_gauge_wave__what_fn,"Gaussian").eq.1) then
H = one - a*dexp(-(x-t)**2/d**2)
end if
C write metric.
if (Minkowski_gauge_wave__diagonal .eq. 1) then
gdxx = half * H + half
gdxy = -half* H + half
gdyy = half * H + half
guxx = half / H + half
guxy =-half / H + half
guyy = half / H + half
else
gdxx = H
gdxy = zero
gdyy = one
guxx = one/H
guxy = zero
guyy = one
end if
gdtt = - H
gdtx = zero
gdty = zero
gdtz = zero
gdzx = zero
gdyz = zero
gdzz = one
C and upper metric.
gutt = - one/H
gutx = zero
guty = zero
gutz = zero
guyz = zero
guzx = zero
guzz = one
return
end
|