testKernel.c 2.41 KB
Newer Older
1
2
/*******************************************************************************
 * This file is part of SWIFT.
3
4
 * Copyright (C) 2016 Matthieu Schaller (matthieu.schaller@durham.ac.uk)
 *                    James Willis (james.s.willis@durham.ac.uk)
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU Lesser General Public License as published
 * by the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 *
 ******************************************************************************/

21
22
23
#define NO__AVX__
#include "vector.h"

24
#include "swift.h"
25
#include "kernel_hydro.h"
26

27
28
#define numPoints 64

29
30
31
int main() {

  const float h = const_eta_kernel;
32
33
  float W[numPoints] = {0.f};
  float dW[numPoints] = {0.f};
34

35
  printf("\nSerial Output\n");
36
  printf("-------------\n");
37

Matthieu Schaller's avatar
Matthieu Schaller committed
38
  for (int i = 0; i < numPoints; ++i) {
39

Matthieu Schaller's avatar
Matthieu Schaller committed
40
    const float x = i * 2.5f / numPoints;
41
    kernel_deval(x / h, &W[i], &dW[i]);
42

Matthieu Schaller's avatar
Matthieu Schaller committed
43
    printf("%2d: h= %f H= %f x=%f W(x,h)=%f dW(x,h)=%f\n", i, h,
44
           h * kernel_gamma, x, W[i], dW[i]);
45
46
  }

Matthieu Schaller's avatar
Matthieu Schaller committed
47
  printf("\nVector Output for VEC_SIZE=%d\n", VEC_SIZE);
48
  printf("-------------\n");
Matthieu Schaller's avatar
Matthieu Schaller committed
49
  for (int i = 0; i < numPoints; i += VEC_SIZE) {
50
51

    vector vx, vx_h;
52
    vector W_vec, dW_vec;
53

Matthieu Schaller's avatar
Matthieu Schaller committed
54
55
    for (int j = 0; j < VEC_SIZE; j++) {
      vx.f[j] = (i + j) * 2.5f / numPoints;
56
    }
Matthieu Schaller's avatar
Matthieu Schaller committed
57

58
59
    vx_h.v = vx.v / vec_set1(h);

60
    kernel_deval_vec(&vx_h, &W_vec, &dW_vec);
61

Matthieu Schaller's avatar
Matthieu Schaller committed
62
63
    for (int j = 0; j < VEC_SIZE; j++) {
      printf("%2d: h= %f H= %f x=%f W(x,h)=%f dW(x,h)=%f\n", i + j, h,
64
65
66
67
68
69
70
71
72
73
74
75
             h * kernel_gamma, vx.f[j], W_vec.f[j], dW_vec.f[j]);

      if (fabsf(W_vec.f[j] - W[i + j]) > 2e-7) {
        printf("Invalid value ! scalar= %e, vector= %e\n", W[i + j],
               W_vec.f[j]);
        return 1;
      }
      if (fabsf(dW_vec.f[j] - dW[i + j]) > 2e-7) {
        printf("Invalid value ! scalar= %e, vector= %e\n", dW[i + j],
               dW_vec.f[j]);
        return 1;
      }
76
77
    }
  }
78
79

  printf("\nAll values are consistent\n");
80
81
  return 0;
}