aboutsummaryrefslogtreecommitdiff
path: root/src/Dissipation_6_5.F90
blob: 34f39f1c33c743a63ec695a260dc9881395d76af (plain)
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
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
! $Header$

#include "cctk.h"

subroutine dissipation_6_5 (var, ni, nj, nk, bb, gsize, delta, epsdis, rhs)

  implicit none

  integer ::  ni, nj, nk
  CCTK_REAL, dimension(ni,nj,nk), intent(in) :: var
  CCTK_REAL, dimension(ni,nj,nk), intent(inout) :: rhs
  CCTK_INT, dimension(6), intent(in) :: bb
  CCTK_INT, dimension(3), intent(in) :: gsize
  CCTK_REAL, dimension(3), intent(in) :: delta
  CCTK_REAL, intent(in) :: epsdis 

  CCTK_REAL :: zero = 0.0
  integer, parameter :: wp = kind(zero)
  integer :: i, j, k
  CCTK_REAL, dimension(4), save :: ac = (/ -20.0_wp, 15.0_wp, -6.0_wp, 1.0_wp /)
  CCTK_REAL, dimension(10,7) :: a
  CCTK_REAL :: idel

  CCTK_INT :: il, ir, jl, jr, kl, kr

!  print*,"dissipation_6_5 called"
!  print*,ni,nj,nk
!  print*,bb
!  print*,gsize
!  print*,delta
!  print*,epsdis
!  stop
  call set_coeff ( a, delta(1) )

  idel = epsdis / ( 64 * delta(1) )

  if ( bb(1) == 0 ) then
    il = 1 + gsize(1)
  else
    rhs(1,:,:) = rhs(1,:,:) + &
                 ( a(1,1) * var(1,:,:) + a(2,1) * var(2,:,:) + &
                   a(3,1) * var(3,:,:) + a(4,1) * var(4,:,:) ) * idel
    rhs(2,:,:) = rhs(2,:,:) + &
                 ( a(1,2) * var(1,:,:) + a(2,2) * var(2,:,:) + &
                   a(3,2) * var(3,:,:) + a(4,2) * var(4,:,:) + &
                   a(5,2) * var(5,:,:) + a(6,2) * var(6,:,:) + &
                   a(7,2) * var(7,:,:) + a(8,2) * var(8,:,:) ) * idel
    rhs(3,:,:) = rhs(3,:,:) + &
                 ( a(1,3) * var(1,:,:) + a(2,3) * var(2,:,:) + &
                   a(3,3) * var(3,:,:) + a(4,3) * var(4,:,:) + &
                   a(5,3) * var(5,:,:) + a(6,3) * var(6,:,:) + &
                   a(7,3) * var(7,:,:) + a(8,3) * var(8,:,:) + &
                   a(9,3) * var(9,:,:) ) * idel
    rhs(4,:,:) = rhs(4,:,:) + &
                 ( a(1,4) * var(1,:,:) + a(2,4) * var(2,:,:) + &
                   a(3,4) * var(3,:,:) + a(4,4) * var(4,:,:) + &
                   a(5,4) * var(5,:,:) + a(6,4) * var(6,:,:) + &
                   a(7,4) * var(7,:,:) + a(8,4) * var(8,:,:) + &
                   a(9,4) * var(9,:,:) + a(10,4) * var(10,:,:) ) * idel
    rhs(5,:,:) = rhs(5,:,:) + &
                 ( a(1,5) * var(1,:,:) + a(2,5) * var(2,:,:) + &
                   a(3,5) * var(3,:,:) + a(4,5) * var(4,:,:) + &
                   a(5,5) * var(5,:,:) + a(6,5) * var(6,:,:) + &
                   a(7,5) * var(7,:,:) + a(8,5) * var(8,:,:) + &
                   a(9,5) * var(9,:,:) + a(10,5) * var(10,:,:) ) * idel
    rhs(6,:,:) = rhs(6,:,:) + &
                 ( a(1,6) * var(1,:,:) + a(2,6) * var(2,:,:) + &
                   a(3,6) * var(3,:,:) + a(4,6) * var(4,:,:) + &
                   a(5,6) * var(5,:,:) + a(6,6) * var(6,:,:) + &
                   a(7,6) * var(7,:,:) + a(8,6) * var(8,:,:) + &
                   a(9,6) * var(9,:,:) + a(10,6) * var(10,:,:) ) * idel
    rhs(7,:,:) = rhs(7,:,:) + &
                 ( a(1,7) * var(1,:,:) + a(2,7) * var(2,:,:) + &
                   a(3,7) * var(3,:,:) + a(4,7) * var(4,:,:) + &
                   a(5,7) * var(5,:,:) + a(6,7) * var(6,:,:) + &
                   a(7,7) * var(7,:,:) + a(8,7) * var(8,:,:) + &
                   a(9,7) * var(9,:,:) + a(10,7) * var(10,:,:) ) * idel
    il = 8
  end if
  if ( bb(2) == 0 ) then
    ir = ni - gsize(1)
  else
    rhs(ni-6,:,:) = rhs(ni-6,:,:) + &
                    ( a(1,7) * var(ni,:,:) + a(2,7) * var(ni-1,:,:) + &
                      a(3,7) * var(ni-2,:,:) + a(4,7) * var(ni-3,:,:) + &
                      a(5,7) * var(ni-4,:,:) + a(6,7) * var(ni-5,:,:) + &
                      a(7,7) * var(ni-6,:,:) + a(8,7) * var(ni-7,:,:) + &
                      a(9,7) * var(ni-8,:,:) + a(10,7) * var(ni-9,:,:) ) * idel
    rhs(ni-5,:,:) = rhs(ni-5,:,:) + &
                    ( a(1,6) * var(ni,:,:) + a(2,6) * var(ni-1,:,:) + &
                      a(3,6) * var(ni-2,:,:) + a(4,6) * var(ni-3,:,:) + &
                      a(5,6) * var(ni-4,:,:) + a(6,6) * var(ni-5,:,:) + &
                      a(7,6) * var(ni-6,:,:) + a(8,6) * var(ni-7,:,:) + &
                      a(9,6) * var(ni-8,:,:) + a(10,6) * var(ni-9,:,:) ) * idel
    rhs(ni-4,:,:) = rhs(ni-4,:,:) + &
                    ( a(1,5) * var(ni,:,:) + a(2,5) * var(ni-1,:,:) + &
                      a(3,5) * var(ni-2,:,:) + a(4,5) * var(ni-3,:,:) + &
                      a(5,5) * var(ni-4,:,:) + a(6,5) * var(ni-5,:,:) + &
                      a(7,5) * var(ni-6,:,:) + a(8,5) * var(ni-7,:,:) + &
                      a(9,5) * var(ni-8,:,:) + a(10,5) * var(ni-9,:,:) ) * idel
    rhs(ni-3,:,:) = rhs(ni-3,:,:) + &
                    ( a(1,4) * var(ni,:,:) + a(2,4) * var(ni-1,:,:) + &
                      a(3,4) * var(ni-2,:,:) + a(4,4) * var(ni-3,:,:) + &
                      a(5,4) * var(ni-4,:,:) + a(6,4) * var(ni-5,:,:) + &
                      a(7,4) * var(ni-6,:,:) + a(8,4) * var(ni-7,:,:) + &
                      a(9,4) * var(ni-8,:,:) + a(10,4) * var(ni-9,:,:) ) * idel
    rhs(ni-2,:,:) = rhs(ni-2,:,:) + &
                    ( a(1,3) * var(ni,:,:) + a(2,3) * var(ni-1,:,:) + &
                      a(3,3) * var(ni-2,:,:) + a(4,3) * var(ni-3,:,:) + &
                      a(5,3) * var(ni-4,:,:) + a(6,3) * var(ni-5,:,:) + &
                      a(7,3) * var(ni-6,:,:) + a(8,3) * var(ni-7,:,:) + &
                      a(9,3) * var(ni-8,:,:) ) * idel
    rhs(ni-1,:,:) = rhs(ni-1,:,:) + &
                    ( a(1,2) * var(ni,:,:) + a(2,2) * var(ni-1,:,:) + &
                      a(3,2) * var(ni-2,:,:) + a(4,2) * var(ni-3,:,:) + &
                      a(5,2) * var(ni-4,:,:) + a(6,2) * var(ni-5,:,:) + &
                      a(7,2) * var(ni-6,:,:) + a(8,2) * var(ni-7,:,:) ) * idel
    rhs(ni,:,:)   = rhs(ni,:,:) + &
                    ( a(1,1) * var(ni,:,:) + a(2,1) * var(ni-1,:,:) + &
                      a(3,1) * var(ni-2,:,:) + a(4,1) * var(ni-3,:,:) ) * idel
    ir = ni - 7
  end if
  rhs(il:ir,:,:) = rhs(il:ir,:,:) + &
                   ( -20.0_wp * var(il:ir,:,:) + &
                      15.0_wp * ( var(il-1:ir-1,:,:) + &
                                  var(il+1:ir+1,:,:) ) - &
                       6.0_wp * ( var(il-2:ir-2,:,:) + &
                                  var(il+2:ir+2,:,:) ) + &
                                ( var(il-3:ir-3,:,:) + &
                                  var(il+3:ir+3,:,:) ) ) * idel

  call set_coeff ( a, delta(2) )

  idel = epsdis / ( 64 * delta(2) )

  if ( bb(3) == 0 ) then
    jl = 1 + gsize(2)
  else
    rhs(:,1,:) = rhs(:,1,:) + &
                 ( a(1,1) * var(:,1,:) + a(2,1) * var(:,2,:) + &
                   a(3,1) * var(:,3,:) + a(4,1) * var(:,4,:) ) * idel
    rhs(:,2,:) = rhs(:,2,:) + &
                 ( a(1,2) * var(:,1,:) + a(2,2) * var(:,2,:) + &
                   a(3,2) * var(:,3,:) + a(4,2) * var(:,4,:) + &
                   a(5,2) * var(:,5,:) + a(6,2) * var(:,6,:) + &
                   a(7,2) * var(:,7,:) + a(8,2) * var(:,8,:) ) * idel
    rhs(:,3,:) = rhs(:,3,:) + &
                 ( a(1,3) * var(:,1,:) + a(2,3) * var(:,2,:) + &
                   a(3,3) * var(:,3,:) + a(4,3) * var(:,4,:) + &
                   a(5,3) * var(:,5,:) + a(6,3) * var(:,6,:) + &
                   a(7,3) * var(:,7,:) + a(8,3) * var(:,8,:) + &
                   a(9,3) * var(:,9,:) ) * idel
    rhs(:,4,:) = rhs(:,4,:) + &
                 ( a(1,4) * var(:,1,:) + a(2,4) * var(:,2,:) + &
                   a(3,4) * var(:,3,:) + a(4,4) * var(:,4,:) + &
                   a(5,4) * var(:,5,:) + a(6,4) * var(:,6,:) + &
                   a(7,4) * var(:,7,:) + a(8,4) * var(:,8,:) + &
                   a(9,4) * var(:,9,:) + a(10,4) * var(:,10,:) ) * idel
    rhs(:,5,:) = rhs(:,5,:) + &
                 ( a(1,5) * var(:,1,:) + a(2,5) * var(:,2,:) + &
                   a(3,5) * var(:,3,:) + a(4,5) * var(:,4,:) + &
                   a(5,5) * var(:,5,:) + a(6,5) * var(:,6,:) + &
                   a(7,5) * var(:,7,:) + a(8,5) * var(:,8,:) + &
                   a(9,5) * var(:,9,:) + a(10,5) * var(:,10,:) ) * idel
    rhs(:,6,:) = rhs(:,6,:) + &
                 ( a(1,6) * var(:,1,:) + a(2,6) * var(:,2,:) + &
                   a(3,6) * var(:,3,:) + a(4,6) * var(:,4,:) + &
                   a(5,6) * var(:,5,:) + a(6,6) * var(:,6,:) + &
                   a(7,6) * var(:,7,:) + a(8,6) * var(:,8,:) + &
                   a(9,6) * var(:,9,:) + a(10,6) * var(:,10,:) ) * idel
    rhs(:,7,:) = rhs(:,7,:) + &
                 ( a(1,7) * var(:,1,:) + a(2,7) * var(:,2,:) + &
                   a(3,7) * var(:,3,:) + a(4,7) * var(:,4,:) + &
                   a(5,7) * var(:,5,:) + a(6,7) * var(:,6,:) + &
                   a(7,7) * var(:,7,:) + a(8,7) * var(:,8,:) + &
                   a(9,7) * var(:,9,:) + a(10,7) * var(:,10,:) ) * idel
    jl = 8
  end if
  if ( bb(4) == 0 ) then
    jr = nj - gsize(2)
  else
    rhs(:,nj-6,:) = rhs(:,nj-6,:) + &
                    ( a(1,7) * var(:,nj,:) + a(2,7) * var(:,nj-1,:) + &
                      a(3,7) * var(:,nj-2,:) + a(4,7) * var(:,nj-3,:) + &
                      a(5,7) * var(:,nj-4,:) + a(6,7) * var(:,nj-5,:) + &
                      a(7,7) * var(:,nj-6,:) + a(8,7) * var(:,nj-7,:) + &
                      a(9,7) * var(:,nj-8,:) + a(10,7) * var(:,nj-9,:) ) * idel
    rhs(:,nj-5,:) = rhs(:,nj-5,:) + &
                    ( a(1,6) * var(:,nj,:) + a(2,6) * var(:,nj-1,:) + &
                      a(3,6) * var(:,nj-2,:) + a(4,6) * var(:,nj-3,:) + &
                      a(5,6) * var(:,nj-4,:) + a(6,6) * var(:,nj-5,:) + &
                      a(7,6) * var(:,nj-6,:) + a(8,6) * var(:,nj-7,:) + &
                      a(9,6) * var(:,nj-8,:) + a(10,6) * var(:,nj-9,:) ) * idel
    rhs(:,nj-4,:) = rhs(:,nj-4,:) + &
                    ( a(1,5) * var(:,nj,:) + a(2,5) * var(:,nj-1,:) + &
                      a(3,5) * var(:,nj-2,:) + a(4,5) * var(:,nj-3,:) + &
                      a(5,5) * var(:,nj-4,:) + a(6,5) * var(:,nj-5,:) + &
                      a(7,5) * var(:,nj-6,:) + a(8,5) * var(:,nj-7,:) + &
                      a(9,5) * var(:,nj-8,:) + a(10,5) * var(:,nj-9,:) ) * idel
    rhs(:,nj-3,:) = rhs(:,nj-3,:) + &
                    ( a(1,4) * var(:,nj,:) + a(2,4) * var(:,nj-1,:) + &
                      a(3,4) * var(:,nj-2,:) + a(4,4) * var(:,nj-3,:) + &
                      a(5,4) * var(:,nj-4,:) + a(6,4) * var(:,nj-5,:) + &
                      a(7,4) * var(:,nj-6,:) + a(8,4) * var(:,nj-7,:) + &
                      a(9,4) * var(:,nj-8,:) + a(10,4) * var(:,nj-9,:) ) * idel
    rhs(:,nj-2,:) = rhs(:,nj-2,:) + &
                    ( a(1,3) * var(:,nj,:) + a(2,3) * var(:,nj-1,:) + &
                      a(3,3) * var(:,nj-2,:) + a(4,3) * var(:,nj-3,:) + &
                      a(5,3) * var(:,nj-4,:) + a(6,3) * var(:,nj-5,:) + &
                      a(7,3) * var(:,nj-6,:) + a(8,3) * var(:,nj-7,:) + &
                      a(9,3) * var(:,nj-8,:) ) * idel
    rhs(:,nj-1,:) = rhs(:,nj-1,:) + &
                    ( a(1,2) * var(:,nj,:) + a(2,2) * var(:,nj-1,:) + &
                      a(3,2) * var(:,nj-2,:) + a(4,2) * var(:,nj-3,:) + &
                      a(5,2) * var(:,nj-4,:) + a(6,2) * var(:,nj-5,:) + &
                      a(7,2) * var(:,nj-6,:) + a(8,2) * var(:,nj-7,:) ) * idel
    rhs(:,nj,:)   = rhs(:,nj,:) + &
                    ( a(1,1) * var(:,nj,:) + a(2,1) * var(:,nj-1,:) + &
                      a(3,1) * var(:,nj-2,:) + a(4,1) * var(:,nj-3,:) ) * idel
    jr = nj - 7
  end if
  rhs(:,jl:jr,:) = rhs(:,jl:jr,:) + &
                   ( -20.0_wp * var(:,jl:jr,:) + &
                      15.0_wp * ( var(:,jl-1:jr-1,:) + &
                                  var(:,jl+1:jr+1,:) ) - &
                       6.0_wp * ( var(:,jl-2:jr-2,:) + &
                                  var(:,jl+2:jr+2,:) ) + &
                                ( var(:,jl-3:jr-3,:) + &
                                  var(:,jl+3:jr+3,:) ) ) * idel

  call set_coeff ( a, delta(3) )

  idel = epsdis / ( 64 * delta(3) )

  if ( bb(5) == 0 ) then
    kl = 1 + gsize(3)
  else
    rhs(:,:,1) = rhs(:,:,1) + &
                 ( a(1,1) * var(:,:,1) + a(2,1) * var(:,:,2) + &
                   a(3,1) * var(:,:,3) + a(4,1) * var(:,:,4) ) * idel
    rhs(:,:,2) = rhs(:,:,2) + &
                 ( a(1,2) * var(:,:,1) + a(2,2) * var(:,:,2) + &
                   a(3,2) * var(:,:,3) + a(4,2) * var(:,:,4) + &
                   a(5,2) * var(:,:,5) + a(6,2) * var(:,:,6) + &
                   a(7,2) * var(:,:,7) + a(8,2) * var(:,:,8) ) * idel
    rhs(:,:,3) = rhs(:,:,3) + &
                 ( a(1,3) * var(:,:,1) + a(2,3) * var(:,:,2) + &
                   a(3,3) * var(:,:,3) + a(4,3) * var(:,:,4) + &
                   a(5,3) * var(:,:,5) + a(6,3) * var(:,:,6) + &
                   a(7,3) * var(:,:,7) + a(8,3) * var(:,:,8) + &
                   a(9,3) * var(:,:,9) ) * idel
    rhs(:,:,4) = rhs(:,:,4) + &
                 ( a(1,4) * var(:,:,1) + a(2,4) * var(:,:,2) + &
                   a(3,4) * var(:,:,3) + a(4,4) * var(:,:,4) + &
                   a(5,4) * var(:,:,5) + a(6,4) * var(:,:,6) + &
                   a(7,4) * var(:,:,7) + a(8,4) * var(:,:,8) + &
                   a(9,4) * var(:,:,9) + a(10,4) * var(:,:,10) ) * idel
    rhs(:,:,5) = rhs(:,:,5) + &
                 ( a(1,5) * var(:,:,1) + a(2,5) * var(:,:,2) + &
                   a(3,5) * var(:,:,3) + a(4,5) * var(:,:,4) + &
                   a(5,5) * var(:,:,5) + a(6,5) * var(:,:,6) + &
                   a(7,5) * var(:,:,7) + a(8,5) * var(:,:,8) + &
                   a(9,5) * var(:,:,9) + a(10,5) * var(:,:,10) ) * idel
    rhs(:,:,6) = rhs(:,:,6) + &
                 ( a(1,6) * var(:,:,1) + a(2,6) * var(:,:,2) + &
                   a(3,6) * var(:,:,3) + a(4,6) * var(:,:,4) + &
                   a(5,6) * var(:,:,5) + a(6,6) * var(:,:,6) + &
                   a(7,6) * var(:,:,7) + a(8,6) * var(:,:,8) + &
                   a(9,6) * var(:,:,9) + a(10,6) * var(:,:,10) ) * idel
    rhs(:,:,7) = rhs(:,:,7) + &
                 ( a(1,7) * var(:,:,1) + a(2,7) * var(:,:,2) + &
                   a(3,7) * var(:,:,3) + a(4,7) * var(:,:,4) + &
                   a(5,7) * var(:,:,5) + a(6,7) * var(:,:,6) + &
                   a(7,7) * var(:,:,7) + a(8,7) * var(:,:,8) + &
                   a(9,7) * var(:,:,9) + a(10,7) * var(:,:,10) ) * idel
    kl = 8
  end if
  if ( bb(6) == 0 ) then
    kr = nk - gsize(3)
  else
    rhs(:,:,nk-6) = rhs(:,:,nk-6) + &
                    ( a(1,7) * var(:,:,nk) + a(2,7) * var(:,:,nk-1) + &
                      a(3,7) * var(:,:,nk-2) + a(4,7) * var(:,:,nk-3) + &
                      a(5,7) * var(:,:,nk-4) + a(6,7) * var(:,:,nk-5) + &
                      a(7,7) * var(:,:,nk-6) + a(8,7) * var(:,:,nk-7) + &
                      a(9,7) * var(:,:,nk-8) + a(10,7) * var(:,:,nk-9) ) * idel
    rhs(:,:,nk-5) = rhs(:,:,nk-5) + &
                    ( a(1,6) * var(:,:,nk) + a(2,6) * var(:,:,nk-1) + &
                      a(3,6) * var(:,:,nk-2) + a(4,6) * var(:,:,nk-3) + &
                      a(5,6) * var(:,:,nk-4) + a(6,6) * var(:,:,nk-5) + &
                      a(7,6) * var(:,:,nk-6) + a(8,6) * var(:,:,nk-7) + &
                      a(9,6) * var(:,:,nk-8) + a(10,6) * var(:,:,nk-9) ) * idel
    rhs(:,:,nk-4) = rhs(:,:,nk-4) + &
                    ( a(1,5) * var(:,:,nk) + a(2,5) * var(:,:,nk-1) + &
                      a(3,5) * var(:,:,nk-2) + a(4,5) * var(:,:,nk-3) + &
                      a(5,5) * var(:,:,nk-4) + a(6,5) * var(:,:,nk-5) + &
                      a(7,5) * var(:,:,nk-6) + a(8,5) * var(:,:,nk-7) + &
                      a(9,5) * var(:,:,nk-8) + a(10,5) * var(:,:,nk-9) ) * idel
    rhs(:,:,nk-3) = rhs(:,:,nk-3) + &
                    ( a(1,4) * var(:,:,nk) + a(2,4) * var(:,:,nk-1) + &
                      a(3,4) * var(:,:,nk-2) + a(4,4) * var(:,:,nk-3) + &
                      a(5,4) * var(:,:,nk-4) + a(6,4) * var(:,:,nk-5) + &
                      a(7,4) * var(:,:,nk-6) + a(8,4) * var(:,:,nk-7) + &
                      a(9,4) * var(:,:,nk-8) + a(10,4) * var(:,:,nk-9) ) * idel
    rhs(:,:,nk-2) = rhs(:,:,nk-2) + &
                    ( a(1,3) * var(:,:,nk) + a(2,3) * var(:,:,nk-1) + &
                      a(3,3) * var(:,:,nk-2) + a(4,3) * var(:,:,nk-3) + &
                      a(5,3) * var(:,:,nk-4) + a(6,3) * var(:,:,nk-5) + &
                      a(7,3) * var(:,:,nk-6) + a(8,3) * var(:,:,nk-7) + &
                      a(9,3) * var(:,:,nk-8) ) * idel
    rhs(:,:,nk-1) = rhs(:,:,nk-1) + &
                    ( a(1,2) * var(:,:,nk) + a(2,2) * var(:,:,nk-1) + &
                      a(3,2) * var(:,:,nk-2) + a(4,2) * var(:,:,nk-3) + &
                      a(5,2) * var(:,:,nk-4) + a(6,2) * var(:,:,nk-5) + &
                      a(7,2) * var(:,:,nk-6) + a(8,2) * var(:,:,nk-7) ) * idel
    rhs(:,:,nk)   = rhs(:,:,nk) + &
                    ( a(1,1) * var(:,:,nk) + a(2,1) * var(:,:,nk-1) + &
                      a(3,1) * var(:,:,nk-2) + a(4,1) * var(:,:,nk-3) ) * idel
    kr = nk - 7
  end if
  rhs(:,:,kl:kr) = rhs(:,:,kl:kr) + &
                   ( -20.0_wp * var(:,:,kl:kr) + &
                      15.0_wp * ( var(:,:,kl-1:kr-1) + &
                                  var(:,:,kl+1:kr+1) ) - &
                       6.0_wp * ( var(:,:,kl-2:kr-2) + &
                                  var(:,:,kl+2:kr+2) ) + &
                                ( var(:,:,kl-3:kr-3) + &
                                  var(:,:,kl+3:kr+3) ) ) * idel
contains

  subroutine set_coeff ( a, dx )

    implicit none

    CCTK_REAL, dimension(10,7), intent(OUT) :: a
    CCTK_REAL, intent(in) :: dx
    CCTK_REAL :: zero = 0.0
    integer, parameter :: wp = kind(zero)
    CCTK_REAL :: b1, b2, b3, b4, b5, b6

    b1 = dx**2; b2 = (2+25*b1)/27; b3 = (7+20*b1)/27;
    b4 = (1+b1)/2; b5 = (20+7*b1)/27; b6 = (25+2*b1)/27;

    a(1,1) = -4.566666666666666666666666666666666666666666666667_wp*b1 &
            - 4.566666666666666666666666666666666666666666666667_wp*b2 &
            - 4.566666666666666666666666666666666666666666666667_wp*b3
    a(2,1) = 13.7_wp*b1 + 13.7_wp*b2 + 13.7_wp*b3
    a(3,1) = -13.7_wp*b1 - 13.7_wp*b2 - 13.7_wp*b3
    a(4,1) = 4.566666666666666666666666666666666666666666666667_wp*b1 &
           + 4.566666666666666666666666666666666666666666666667_wp*b2 &
           + 4.566666666666666666666666666666666666666666666667_wp*b3
    a(5,1) = zero
    a(6,1) = zero
    a(7,1) = zero
    a(8,1) = zero
    a(9,1) = zero
    a(10,1) = zero
    
    a(1,2) = 0.797987233606384804466079266492184563887676357853_wp*b1 &
           + 0.797987233606384804466079266492184563887676357853_wp*b2 &
           + 0.797987233606384804466079266492184563887676357853_wp*b3
    a(2,2) = -2.393961700819154413398237799476553691663029073559_wp*b1 &
            - 2.393961700819154413398237799476553691663029073559_wp*b2 &
            - 2.393961700819154413398237799476553691663029073559_wp*b3 &
            + 0.307359194675104101198317748428244417750742420031_wp*b4
    a(3,2) = 2.393961700819154413398237799476553691663029073559_wp*b1 &
           + 2.393961700819154413398237799476553691663029073559_wp*b2 &
           + 2.393961700819154413398237799476553691663029073559_wp*b3 &
           - 0.922077584025312303594953245284733253252227260092_wp*b4 &
           - 0.0336105345682464004522484486235834134456512678_wp*b5
    a(4,2) = -0.797987233606384804466079266492184563887676357853_wp*b1 &
            - 0.797987233606384804466079266492184563887676357853_wp*b2 &
            - 0.797987233606384804466079266492184563887676357853_wp*b3 &
            + 0.922077584025312303594953245284733253252227260092_wp*b4 &
            + 0.1008316037047392013567453458707502403369538034_wp*b5 &
            - 0.2929089778363370916823203340749080588489943555_wp*b6
    a(5,2) = 0.09404587787607681409837743471424421680278209869_wp &
           - 0.307359194675104101198317748428244417750742420031_wp*b4 &
           - 0.1008316037047392013567453458707502403369538034_wp*b5 &
           + 0.878726933509011275046961002224724176546983066499_wp*b6
    a(6,2) = -0.282137633628230442295132304142732650409647604264_wp &
            + 0.0336105345682464004522484486235834134456512678_wp*b5 &
            - 0.878726933509011275046961002224724176546983066499_wp*b6
    a(7,2) = 0.282137633628230442295132304142732650412591061208_wp &
           + 0.2929089778363370916823203340749080588489943555_wp*b6
    a(8,2) = -0.094045877876076814098377434714244216807708544941_wp
    a(9,2) = zero
    a(10,2) = zero
    
    a(1,3) = -0.570050756575490634640420716139511539225905452034_wp*b1 &
            - 0.570050756575490634640420716139511539225905452034_wp*b2 &
            - 0.570050756575490634640420716139511539225905452034_wp*b3
    a(2,3) = 1.710152269726471903921262148418534617677716356102_wp*b1 &
           + 1.710152269726471903921262148418534617677716356102_wp*b2 &
           + 1.710152269726471903921262148418534617677716356102_wp*b3 &
           + 0.65012922388697115376186819507780973585422870837_wp*b4
    a(3,3) = -1.710152269726471903921262148418534617677716356102_wp*b1 &
            - 1.710152269726471903921262148418534617677716356102_wp*b2 &
            - 1.710152269726471903921262148418534617677716356102_wp*b3 &
            - 1.950387671660913461285604585233429207562686125111_wp*b4 &
            + 0.222968626602313512338688373037925009289627989014_wp*b5
    a(4,3) = 0.570050756575490634640420716139511539225905452034_wp*b1 &
           + 0.570050756575490634640420716139511539225905452034_wp*b2 &
           + 0.570050756575490634640420716139511539225905452034_wp*b3 &
           + 1.950387671660913461285604585233429207562686125111_wp*b4 &
           - 0.668905879806940537016065119113775027868883967041_wp*b5 &
           - 0.10042075802449997038657833546993392801822845339_wp*b6
    a(5,3) = -0.223539235539512073697159351113375565174639048593_wp &
            - 0.65012922388697115376186819507780973585422870837_wp*b4 &
            + 0.668905879806940537016065119113775027868883967041_wp*b5 &
            + 0.30126227407349991115973500640980178405468536017_wp*b6
    a(6,3) = 0.74867234655369702626509041115338883998616371461_wp &
           - 0.222968626602313512338688373037925009289627989014_wp*b5 &
           - 0.30126227407349991115973500640980178405468536017_wp*b6
    a(7,3) = -0.904781626424018636612315126779913128909937434049_wp &
            + 0.10042075802449997038657833546993392801822845339_wp*b6
    a(8,3) = 0.45770315534499448921799642455316199855922050042_wp
    a(9,3) = -0.078054639935160805173612357813262144460088314167_wp
    a(10,3) = zero
    
    a(1,4) = -0.432850298736358291696660164313057377022426239952_wp*b1 &
            - 0.432850298736358291696660164313057377022426239952_wp*b2 &
            - 0.432850298736358291696660164313057377022426239952_wp*b3
    a(2,4) = 1.298550896209074875089980492939172131067278719856_wp*b1 &
           + 1.298550896209074875089980492939172131067278719856_wp*b2 &
           + 1.298550896209074875089980492939172131067278719856_wp*b3 &
           - 0.13602180681849030343231048896787005715681064672_wp*b4
    a(3,4) = -1.298550896209074875089980492939172131067278719856_wp*b1 &
            - 1.298550896209074875089980492939172131067278719856_wp*b2 &
            - 1.298550896209074875089980492939172131067278719856_wp*b3 &
            + 0.40806542045547091029693146690361017147043194016_wp*b4 &
            + 0.139708407134619223127669197053218276788706051412_wp*b5
    a(4,4) = 0.432850298736358291696660164313057377022426239952_wp*b1 &
           + 0.432850298736358291696660164313057377022426239952_wp*b2 &
           + 0.432850298736358291696660164313057377022426239952_wp*b3 &
           - 0.40806542045547091029693146690361017147043194016_wp*b4 &
           - 0.419125221403857669383007591159654830366118154235_wp*b5 &
           + 0.422237850757834010213265187828443211383102779457_wp*b6
    a(5,4) = -0.231917524240705800245249424344181664391294057734_wp &
            + 0.13602180681849030343231048896787005715681064672_wp*b4 &
            + 0.419125221403857669383007591159654830366118154235_wp*b5 &
            - 1.26671355227350203063979556348532963414930833837_wp*b6
    a(6,4) = 0.632304366383038264396056896345266141026261929548_wp &
           - 0.139708407134619223127669197053218276788706051412_wp*b5 &
           + 1.26671355227350203063979556348532963414930833837_wp*b6
    a(7,4) = -0.484284312078557755005703255606118914056056082961_wp &
            - 0.422237850757834010213265187828443211383102779457_wp*b6
    a(8,4) = -0.021798019655498318906737367811423460076462751066_wp
    a(9,4) = 0.126819131218045846472604038781047420172516321492_wp
    a(10,4) = -0.021123641626322236710970887364589522674965359279_wp
    
    a(1,5) = 0.04525671266969037217439909305498103903084448236_wp*b1 &
           + 0.04525671266969037217439909305498103903084448236_wp*b2 &
           + 0.04525671266969037217439909305498103903084448236_wp*b3
    a(2,5) = -0.135770138009071116523197279164943117092533447081_wp*b1 &
            - 0.135770138009071116523197279164943117092533447081_wp*b2 &
            - 0.135770138009071116523197279164943117092533447081_wp*b3 &
            - 0.356074970503407919432108002084434090165232256258_wp*b4
    a(3,5) = 0.135770138009071116523197279164943117092533447081_wp*b1 &
           + 0.135770138009071116523197279164943117092533447081_wp*b2 &
           + 0.135770138009071116523197279164943117092533447081_wp*b3 &
           + 1.068224911510223758296324006253302270495696768774_wp*b4 &
           - 0.161252732964068670989303604260577833997844116999_wp*b5
    a(4,5) = -0.04525671266969037217439909305498103903084448236_wp*b1 &
            - 0.04525671266969037217439909305498103903084448236_wp*b2 &
            - 0.04525671266969037217439909305498103903084448236_wp*b3 &
            - 1.068224911510223758296324006253302270495696768774_wp*b4 &
            + 0.483758198892206012967910812781733501993532350997_wp*b5 &
            + 0.000913145598627662444167435877738793583486988286_wp*b6
    a(5,5) = 0.778188619899941335162434896805956583938706965447_wp &
           + 0.356074970503407919432108002084434090165232256258_wp*b4 &
           - 0.483758198892206012967910812781733501993532350997_wp*b5 &
           - 0.002739436795882987332502307633216380750460964857_wp*b6
    a(6,5) = -2.888060650467623805586479517528702488708606593183_wp &
            + 0.161252732964068670989303604260577833997844116999_wp*b5 &
            + 0.002739436795882987332502307633216380750460964857_wp*b6
    a(7,5) = 4.064047176584057140151337032365458641610672353212_wp &
           - 0.000913145598627662444167435877738793583486988286_wp*b6
    a(8,5) = -2.645663825945841938559482959983726831967447155003_wp
    a(9,5) = 0.760485624510301003198698408956104774243768795869_wp
    a(10,5) = -0.068996944580833734366507860615090679117094366342_wp
    
    a(1,6) = 0.234241201265594841727615788033296717516774995874_wp*b1 &
           + 0.234241201265594841727615788033296717516774995874_wp*b2 &
           + 0.234241201265594841727615788033296717516774995874_wp*b3
    a(2,6) = -0.702723603796784525182847364099890152550324987621_wp*b1 &
            - 0.702723603796784525182847364099890152550324987621_wp*b2 &
            - 0.702723603796784525182847364099890152550324987621_wp*b3 &
            + 0.112108192839478902858478028045243413779536843362_wp*b4
    a(3,6) = 0.702723603796784525182847364099890152550324987621_wp*b1 &
           + 0.702723603796784525182847364099890152550324987621_wp*b2 &
           + 0.702723603796784525182847364099890152550324987621_wp*b3 &
           - 0.336324578518436708575434084135730241338610530087_wp*b4 &
           - 0.311645621737838510439662177349969575726926106823_wp*b5
    a(4,6) = -0.234241201265594841727615788033296717516774995874_wp*b1 &
            - 0.234241201265594841727615788033296717516774995874_wp*b2 &
            - 0.234241201265594841727615788033296717516774995874_wp*b3 &
            + 0.336324578518436708575434084135730241338610530087_wp*b4 &
            + 0.934936865213515531318986532049908727180778320468_wp*b5 &
            - 0.814776445474468482860158311580599018777969910384_wp*b6
    a(5,6) = 1.318650558110101765365290411968234643910731822102_wp &
           - 0.112108192839478902858478028045243413779536843362_wp*b4 &
           - 0.934936865213515531318986532049908727180778320468_wp*b5 &
           + 2.444329336423405448580474934741797056333909731151_wp*b6
    a(6,6) = -4.321946659420587264930126779432660223055921154843_wp &
            + 0.311645621737838510439662177349969575726926106823_wp*b5 &
            - 2.444329336423405448580474934741797056333909731151_wp*b6
    a(7,6) = 4.959408369623725459954942028713585060789912972866_wp &
           + 0.814776445474468482860158311580599018777969910384_wp*b6
    a(8,6) = -2.133050733448670443936969529227140283141530210559_wp
    a(9,6) = 0.082410205158004740903168030202993056583347011382_wp
    a(10,6) = 0.094528259977425742643695837774987744913459559052_wp
    
    a(1,7) = -0.021123641626322236710970887364589522673829626284_wp*b1 &
            - 0.021123641626322236710970887364589522673829626284_wp*b2 &
            - 0.021123641626322236710970887364589522673829626284_wp*b3
    
    a(2,7) = 0.063370924878966710132912662093768568021488878853_wp*b1 &
           + 0.063370924878966710132912662093768568021488878853_wp*b2 &
           + 0.063370924878966710132912662093768568021488878853_wp*b3 &
           - 0.005626019701867024233595198521322111094015702611_wp*b4
    a(3,7) = -0.063370924878966710132912662093768568021488878853_wp*b1 &
            - 0.063370924878966710132912662093768568021488878853_wp*b2 &
            - 0.063370924878966710132912662093768568021488878853_wp*b3 &
            + 0.016878059105601072700785595563966333282047107834_wp*b4 &
            + 0.23814816884096023561030675752649121424056599846_wp*b5
    a(4,7) = 0.021123641626322236710970887364589522673829626284_wp*b1 &
           + 0.021123641626322236710970887364589522673829626284_wp*b2 &
           + 0.021123641626322236710970887364589522673829626284_wp*b3 &
           - 0.016878059105601072700785595563966333282047107834_wp*b4 &
           - 0.71444450652288070683092027257947364272169799538_wp*b5 &
           + 0.504229942386996567949909229607651725484488940042_wp*b6
    a(5,7) = -2.568464018793247324511052938299838510846083027331_wp &
            + 0.005626019701867024233595198521322111094015702611_wp*b4 &
            + 0.71444450652288070683092027257947364272169799538_wp*b5 &
            - 1.512689827160989703849727688822955176453466820125_wp*b6
    a(6,7) = 10.531909539708674517698111289364420212328345593768_wp &
            - 0.23814816884096023561030675752649121424056599846_wp*b5 &
            + 1.512689827160989703849727688822955176453466820125_wp*b6
    a(7,7) = -17.15862642080199236829756567570752704681_wp &
             - 0.504229942386996567949909229607651725484488940042_wp*b6
    a(8,7) = 13.96906221208640324381455867393444497492_wp
    a(9,7) = -5.747563226635290830973600786704797104493652582603_wp
    a(10,7) = 0.973681914435452762269549437413297474901185356942_wp
  end subroutine set_coeff

end subroutine dissipation_6_5