xref: /petsc/src/mat/tests/ex66.c (revision a69119a591a03a9d906b29c0a4e9802e4d7c9795)
1 static char help[] = "Tests MATH2OPUS\n\n";
2 
3 #include <petscmat.h>
4 #include <petscsf.h>
5 
6 static PetscScalar GenEntry_Symm(PetscInt sdim, PetscReal x[], PetscReal y[], void *ctx) {
7   PetscInt  d;
8   PetscReal clength = sdim == 3 ? 0.2 : 0.1;
9   PetscReal dist, diff = 0.0;
10 
11   for (d = 0; d < sdim; d++) diff += (x[d] - y[d]) * (x[d] - y[d]);
12   dist = PetscSqrtReal(diff);
13   return PetscExpReal(-dist / clength);
14 }
15 
16 static PetscScalar GenEntry_Unsymm(PetscInt sdim, PetscReal x[], PetscReal y[], void *ctx) {
17   PetscInt  d;
18   PetscReal clength = sdim == 3 ? 0.2 : 0.1;
19   PetscReal dist, diff = 0.0, nx = 0.0, ny = 0.0;
20 
21   for (d = 0; d < sdim; d++) nx += x[d] * x[d];
22   for (d = 0; d < sdim; d++) ny += y[d] * y[d];
23   for (d = 0; d < sdim; d++) diff += (x[d] - y[d]) * (x[d] - y[d]);
24   dist = PetscSqrtReal(diff);
25   return nx > ny ? PetscExpReal(-dist / clength) : PetscExpReal(-dist / clength) + 1.;
26 }
27 
28 int main(int argc, char **argv) {
29   Mat          A, B, C, D;
30   Vec          v, x, y, Ax, Ay, Bx, By;
31   PetscRandom  r;
32   PetscLayout  map;
33   PetscScalar *Adata = NULL, *Cdata = NULL, scale = 1.0;
34   PetscReal   *coords, nA, nD, nB, err, nX, norms[3];
35   PetscInt     N, n = 64, dim = 1, i, j, nrhs = 11, lda = 0, ldc = 0, ldu = 0, nlr = 7, nt, ntrials = 2;
36   PetscMPIInt  size, rank;
37   PetscBool    testlayout = PETSC_FALSE, flg, symm = PETSC_FALSE, Asymm = PETSC_TRUE, kernel = PETSC_TRUE;
38   PetscBool    checkexpl = PETSC_FALSE, agpu = PETSC_FALSE, bgpu = PETSC_FALSE, cgpu = PETSC_FALSE, flgglob;
39   PetscBool    testtrans, testnorm, randommat = PETSC_TRUE, testorthog, testcompress, testhlru;
40   void (*approxnormfunc)(void);
41   void (*Anormfunc)(void);
42 
43 #if defined(PETSC_HAVE_MPI_INIT_THREAD)
44   PETSC_MPI_THREAD_REQUIRED = MPI_THREAD_MULTIPLE;
45 #endif
46   PetscFunctionBeginUser;
47   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
48   PetscCall(PetscOptionsGetInt(NULL, NULL, "-ng", &N, &flgglob));
49   PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
50   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nrhs", &nrhs, NULL));
51   PetscCall(PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL));
52   PetscCall(PetscOptionsGetInt(NULL, NULL, "-lda", &lda, NULL));
53   PetscCall(PetscOptionsGetInt(NULL, NULL, "-ldc", &ldc, NULL));
54   PetscCall(PetscOptionsGetInt(NULL, NULL, "-nlr", &nlr, NULL));
55   PetscCall(PetscOptionsGetInt(NULL, NULL, "-ldu", &ldu, NULL));
56   PetscCall(PetscOptionsGetInt(NULL, NULL, "-matmattrials", &ntrials, NULL));
57   PetscCall(PetscOptionsGetBool(NULL, NULL, "-randommat", &randommat, NULL));
58   if (!flgglob) PetscCall(PetscOptionsGetBool(NULL, NULL, "-testlayout", &testlayout, NULL));
59   PetscCall(PetscOptionsGetBool(NULL, NULL, "-Asymm", &Asymm, NULL));
60   PetscCall(PetscOptionsGetBool(NULL, NULL, "-symm", &symm, NULL));
61   PetscCall(PetscOptionsGetBool(NULL, NULL, "-kernel", &kernel, NULL));
62   PetscCall(PetscOptionsGetBool(NULL, NULL, "-checkexpl", &checkexpl, NULL));
63   PetscCall(PetscOptionsGetBool(NULL, NULL, "-agpu", &agpu, NULL));
64   PetscCall(PetscOptionsGetBool(NULL, NULL, "-bgpu", &bgpu, NULL));
65   PetscCall(PetscOptionsGetBool(NULL, NULL, "-cgpu", &cgpu, NULL));
66   PetscCall(PetscOptionsGetScalar(NULL, NULL, "-scale", &scale, NULL));
67   if (!Asymm) symm = PETSC_FALSE;
68 
69   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
70 
71   /* Disable tests for unimplemented variants */
72   testtrans    = (PetscBool)(size == 1 || symm);
73   testnorm     = (PetscBool)(size == 1 || symm);
74   testorthog   = (PetscBool)(size == 1 || symm);
75   testcompress = (PetscBool)(size == 1 || symm);
76   testhlru     = (PetscBool)(size == 1);
77 
78   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
79   PetscCall(PetscLayoutCreate(PETSC_COMM_WORLD, &map));
80   if (testlayout) {
81     if (rank % 2) n = PetscMax(2 * n - 5 * rank, 0);
82     else n = 2 * n + rank;
83   }
84   if (!flgglob) {
85     PetscCall(PetscLayoutSetLocalSize(map, n));
86     PetscCall(PetscLayoutSetUp(map));
87     PetscCall(PetscLayoutGetSize(map, &N));
88   } else {
89     PetscCall(PetscLayoutSetSize(map, N));
90     PetscCall(PetscLayoutSetUp(map));
91     PetscCall(PetscLayoutGetLocalSize(map, &n));
92   }
93   PetscCall(PetscLayoutDestroy(&map));
94 
95   if (lda) PetscCall(PetscMalloc1(N * (n + lda), &Adata));
96   PetscCall(MatCreateDense(PETSC_COMM_WORLD, n, n, N, N, Adata, &A));
97   PetscCall(MatDenseSetLDA(A, n + lda));
98 
99   /* Create random points; these are replicated in order to populate a dense matrix and to compare sequential and dense runs
100      The constructor for MATH2OPUS can take as input the distributed coordinates and replicates them internally in case
101      a kernel construction is requested */
102   PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &r));
103   PetscCall(PetscRandomSetFromOptions(r));
104   PetscCall(PetscRandomSetSeed(r, 123456));
105   PetscCall(PetscRandomSeed(r));
106   PetscCall(PetscMalloc1(N * dim, &coords));
107   PetscCall(PetscRandomGetValuesReal(r, N * dim, coords));
108   PetscCall(PetscRandomDestroy(&r));
109 
110   if (kernel || !randommat) {
111     MatH2OpusKernel k = Asymm ? GenEntry_Symm : GenEntry_Unsymm;
112     PetscInt        ist, ien;
113 
114     PetscCall(MatGetOwnershipRange(A, &ist, &ien));
115     for (i = ist; i < ien; i++) {
116       for (j = 0; j < N; j++) PetscCall(MatSetValue(A, i, j, (*k)(dim, coords + i * dim, coords + j * dim, NULL), INSERT_VALUES));
117     }
118     PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
119     PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
120     if (kernel) {
121       PetscCall(MatCreateH2OpusFromKernel(PETSC_COMM_WORLD, n, n, N, N, dim, coords + ist * dim, PETSC_TRUE, k, NULL, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, &B));
122     } else {
123       PetscCall(MatCreateH2OpusFromMat(A, dim, coords + ist * dim, PETSC_TRUE, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, &B));
124     }
125   } else {
126     PetscInt ist;
127 
128     PetscCall(MatGetOwnershipRange(A, &ist, NULL));
129     PetscCall(MatSetRandom(A, NULL));
130     if (Asymm) {
131       PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &B));
132       PetscCall(MatAXPY(A, 1.0, B, SAME_NONZERO_PATTERN));
133       PetscCall(MatDestroy(&B));
134       PetscCall(MatSetOption(A, MAT_SYMMETRIC, PETSC_TRUE));
135     }
136     PetscCall(MatCreateH2OpusFromMat(A, dim, coords + ist * dim, PETSC_TRUE, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, &B));
137   }
138   PetscCall(PetscFree(coords));
139   if (agpu) PetscCall(MatConvert(A, MATDENSECUDA, MAT_INPLACE_MATRIX, &A));
140   PetscCall(MatViewFromOptions(A, NULL, "-A_view"));
141 
142   PetscCall(MatSetOption(B, MAT_SYMMETRIC, symm));
143 
144   /* assemble the H-matrix */
145   PetscCall(MatBindToCPU(B, (PetscBool)!bgpu));
146   PetscCall(MatSetFromOptions(B));
147   PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
148   PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
149   PetscCall(MatViewFromOptions(B, NULL, "-B_view"));
150 
151   /* Test MatScale */
152   PetscCall(MatScale(A, scale));
153   PetscCall(MatScale(B, scale));
154 
155   /* Test MatMult */
156   PetscCall(MatCreateVecs(A, &Ax, &Ay));
157   PetscCall(MatCreateVecs(B, &Bx, &By));
158   PetscCall(VecSetRandom(Ax, NULL));
159   PetscCall(VecCopy(Ax, Bx));
160   PetscCall(MatMult(A, Ax, Ay));
161   PetscCall(MatMult(B, Bx, By));
162   PetscCall(VecViewFromOptions(Ay, NULL, "-mult_vec_view"));
163   PetscCall(VecViewFromOptions(By, NULL, "-mult_vec_view"));
164   PetscCall(VecNorm(Ay, NORM_INFINITY, &nX));
165   PetscCall(VecAXPY(Ay, -1.0, By));
166   PetscCall(VecViewFromOptions(Ay, NULL, "-mult_vec_view"));
167   PetscCall(VecNorm(Ay, NORM_INFINITY, &err));
168   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMult err %g\n", err / nX));
169   PetscCall(VecScale(By, -1.0));
170   PetscCall(MatMultAdd(B, Bx, By, By));
171   PetscCall(VecNorm(By, NORM_INFINITY, &err));
172   PetscCall(VecViewFromOptions(By, NULL, "-mult_vec_view"));
173   if (err > 10. * PETSC_SMALL) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMultAdd err %g\n", err));
174 
175   /* Test MatNorm */
176   PetscCall(MatNorm(A, NORM_INFINITY, &norms[0]));
177   PetscCall(MatNorm(A, NORM_1, &norms[1]));
178   norms[2] = -1.; /* NORM_2 not supported */
179   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "A Matrix norms:        infty=%g, norm_1=%g, norm_2=%g\n", (double)norms[0], (double)norms[1], (double)norms[2]));
180   PetscCall(MatGetOperation(A, MATOP_NORM, &Anormfunc));
181   PetscCall(MatGetOperation(B, MATOP_NORM, &approxnormfunc));
182   PetscCall(MatSetOperation(A, MATOP_NORM, approxnormfunc));
183   PetscCall(MatNorm(A, NORM_INFINITY, &norms[0]));
184   PetscCall(MatNorm(A, NORM_1, &norms[1]));
185   PetscCall(MatNorm(A, NORM_2, &norms[2]));
186   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "A Approx Matrix norms: infty=%g, norm_1=%g, norm_2=%g\n", (double)norms[0], (double)norms[1], (double)norms[2]));
187   if (testnorm) {
188     PetscCall(MatNorm(B, NORM_INFINITY, &norms[0]));
189     PetscCall(MatNorm(B, NORM_1, &norms[1]));
190     PetscCall(MatNorm(B, NORM_2, &norms[2]));
191   } else {
192     norms[0] = -1.;
193     norms[1] = -1.;
194     norms[2] = -1.;
195   }
196   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "B Approx Matrix norms: infty=%g, norm_1=%g, norm_2=%g\n", (double)norms[0], (double)norms[1], (double)norms[2]));
197   PetscCall(MatSetOperation(A, MATOP_NORM, Anormfunc));
198 
199   /* Test MatDuplicate */
200   PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &D));
201   PetscCall(MatSetOption(D, MAT_SYMMETRIC, symm));
202   PetscCall(MatMultEqual(B, D, 10, &flg));
203   if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMult error after MatDuplicate\n"));
204   if (testtrans) {
205     PetscCall(MatMultTransposeEqual(B, D, 10, &flg));
206     if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMultTranspose error after MatDuplicate\n"));
207   }
208   PetscCall(MatDestroy(&D));
209 
210   if (testtrans) { /* MatMultTranspose for nonsymmetric matrices not implemented */
211     PetscCall(VecSetRandom(Ay, NULL));
212     PetscCall(VecCopy(Ay, By));
213     PetscCall(MatMultTranspose(A, Ay, Ax));
214     PetscCall(MatMultTranspose(B, By, Bx));
215     PetscCall(VecViewFromOptions(Ax, NULL, "-multtrans_vec_view"));
216     PetscCall(VecViewFromOptions(Bx, NULL, "-multtrans_vec_view"));
217     PetscCall(VecNorm(Ax, NORM_INFINITY, &nX));
218     PetscCall(VecAXPY(Ax, -1.0, Bx));
219     PetscCall(VecViewFromOptions(Ax, NULL, "-multtrans_vec_view"));
220     PetscCall(VecNorm(Ax, NORM_INFINITY, &err));
221     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMultTranspose err %g\n", err / nX));
222     PetscCall(VecScale(Bx, -1.0));
223     PetscCall(MatMultTransposeAdd(B, By, Bx, Bx));
224     PetscCall(VecNorm(Bx, NORM_INFINITY, &err));
225     if (err > 10. * PETSC_SMALL) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMultTransposeAdd err %g\n", err));
226   }
227   PetscCall(VecDestroy(&Ax));
228   PetscCall(VecDestroy(&Ay));
229   PetscCall(VecDestroy(&Bx));
230   PetscCall(VecDestroy(&By));
231 
232   /* Test MatMatMult */
233   if (ldc) PetscCall(PetscMalloc1(nrhs * (n + ldc), &Cdata));
234   PetscCall(MatCreateDense(PETSC_COMM_WORLD, n, PETSC_DECIDE, N, nrhs, Cdata, &C));
235   PetscCall(MatDenseSetLDA(C, n + ldc));
236   PetscCall(MatSetRandom(C, NULL));
237   if (cgpu) PetscCall(MatConvert(C, MATDENSECUDA, MAT_INPLACE_MATRIX, &C));
238   for (nt = 0; nt < ntrials; nt++) {
239     PetscCall(MatMatMult(B, C, nt ? MAT_REUSE_MATRIX : MAT_INITIAL_MATRIX, PETSC_DEFAULT, &D));
240     PetscCall(MatViewFromOptions(D, NULL, "-bc_view"));
241     PetscCall(PetscObjectBaseTypeCompareAny((PetscObject)D, &flg, MATSEQDENSE, MATMPIDENSE, ""));
242     if (flg) {
243       PetscCall(MatCreateVecs(B, &x, &y));
244       PetscCall(MatCreateVecs(D, NULL, &v));
245       for (i = 0; i < nrhs; i++) {
246         PetscCall(MatGetColumnVector(D, v, i));
247         PetscCall(MatGetColumnVector(C, x, i));
248         PetscCall(MatMult(B, x, y));
249         PetscCall(VecAXPY(y, -1.0, v));
250         PetscCall(VecNorm(y, NORM_INFINITY, &err));
251         if (err > 10. * PETSC_SMALL) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMat err %" PetscInt_FMT " %g\n", i, err));
252       }
253       PetscCall(VecDestroy(&y));
254       PetscCall(VecDestroy(&x));
255       PetscCall(VecDestroy(&v));
256     }
257   }
258   PetscCall(MatDestroy(&D));
259 
260   /* Test MatTransposeMatMult */
261   if (testtrans) { /* MatMultTranspose for nonsymmetric matrices not implemented */
262     for (nt = 0; nt < ntrials; nt++) {
263       PetscCall(MatTransposeMatMult(B, C, nt ? MAT_REUSE_MATRIX : MAT_INITIAL_MATRIX, PETSC_DEFAULT, &D));
264       PetscCall(MatViewFromOptions(D, NULL, "-btc_view"));
265       PetscCall(PetscObjectBaseTypeCompareAny((PetscObject)D, &flg, MATSEQDENSE, MATMPIDENSE, ""));
266       if (flg) {
267         PetscCall(MatCreateVecs(B, &y, &x));
268         PetscCall(MatCreateVecs(D, NULL, &v));
269         for (i = 0; i < nrhs; i++) {
270           PetscCall(MatGetColumnVector(D, v, i));
271           PetscCall(MatGetColumnVector(C, x, i));
272           PetscCall(MatMultTranspose(B, x, y));
273           PetscCall(VecAXPY(y, -1.0, v));
274           PetscCall(VecNorm(y, NORM_INFINITY, &err));
275           if (err > 10. * PETSC_SMALL) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatTransMat err %" PetscInt_FMT " %g\n", i, err));
276         }
277         PetscCall(VecDestroy(&y));
278         PetscCall(VecDestroy(&x));
279         PetscCall(VecDestroy(&v));
280       }
281     }
282     PetscCall(MatDestroy(&D));
283   }
284 
285   /* Test basis orthogonalization */
286   if (testorthog) {
287     PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &D));
288     PetscCall(MatSetOption(D, MAT_SYMMETRIC, symm));
289     PetscCall(MatH2OpusOrthogonalize(D));
290     PetscCall(MatMultEqual(B, D, 10, &flg));
291     if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMult error after basis ortogonalization\n"));
292     PetscCall(MatDestroy(&D));
293   }
294 
295   /* Test matrix compression */
296   if (testcompress) {
297     PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &D));
298     PetscCall(MatSetOption(D, MAT_SYMMETRIC, symm));
299     PetscCall(MatH2OpusCompress(D, PETSC_SMALL));
300     PetscCall(MatDestroy(&D));
301   }
302 
303   /* Test low-rank update */
304   if (testhlru) {
305     Mat          U, V;
306     PetscScalar *Udata = NULL, *Vdata = NULL;
307 
308     if (ldu) {
309       PetscCall(PetscMalloc1(nlr * (n + ldu), &Udata));
310       PetscCall(PetscMalloc1(nlr * (n + ldu + 2), &Vdata));
311     }
312     PetscCall(MatDuplicate(B, MAT_COPY_VALUES, &D));
313     PetscCall(MatCreateDense(PetscObjectComm((PetscObject)D), n, PETSC_DECIDE, N, nlr, Udata, &U));
314     PetscCall(MatDenseSetLDA(U, n + ldu));
315     PetscCall(MatCreateDense(PetscObjectComm((PetscObject)D), n, PETSC_DECIDE, N, nlr, Vdata, &V));
316     if (ldu) PetscCall(MatDenseSetLDA(V, n + ldu + 2));
317     PetscCall(MatSetRandom(U, NULL));
318     PetscCall(MatSetRandom(V, NULL));
319     PetscCall(MatH2OpusLowRankUpdate(D, U, V, 0.5));
320     PetscCall(MatH2OpusLowRankUpdate(D, U, V, -0.5));
321     PetscCall(MatMultEqual(B, D, 10, &flg));
322     if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatMult error after low-rank update\n"));
323     PetscCall(MatDestroy(&D));
324     PetscCall(MatDestroy(&U));
325     PetscCall(PetscFree(Udata));
326     PetscCall(MatDestroy(&V));
327     PetscCall(PetscFree(Vdata));
328   }
329 
330   /* check explicit operator */
331   if (checkexpl) {
332     Mat Be, Bet;
333 
334     PetscCall(MatComputeOperator(B, MATDENSE, &D));
335     PetscCall(MatDuplicate(D, MAT_COPY_VALUES, &Be));
336     PetscCall(MatNorm(D, NORM_FROBENIUS, &nB));
337     PetscCall(MatViewFromOptions(D, NULL, "-expl_view"));
338     PetscCall(MatAXPY(D, -1.0, A, SAME_NONZERO_PATTERN));
339     PetscCall(MatViewFromOptions(D, NULL, "-diff_view"));
340     PetscCall(MatNorm(D, NORM_FROBENIUS, &nD));
341     PetscCall(MatNorm(A, NORM_FROBENIUS, &nA));
342     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Approximation error %g (%g / %g, %g)\n", nD / nA, nD, nA, nB));
343     PetscCall(MatDestroy(&D));
344 
345     if (testtrans) { /* MatMultTranspose for nonsymmetric matrices not implemented */
346       PetscCall(MatTranspose(A, MAT_INPLACE_MATRIX, &A));
347       PetscCall(MatComputeOperatorTranspose(B, MATDENSE, &D));
348       PetscCall(MatDuplicate(D, MAT_COPY_VALUES, &Bet));
349       PetscCall(MatNorm(D, NORM_FROBENIUS, &nB));
350       PetscCall(MatViewFromOptions(D, NULL, "-expl_trans_view"));
351       PetscCall(MatAXPY(D, -1.0, A, SAME_NONZERO_PATTERN));
352       PetscCall(MatViewFromOptions(D, NULL, "-diff_trans_view"));
353       PetscCall(MatNorm(D, NORM_FROBENIUS, &nD));
354       PetscCall(MatNorm(A, NORM_FROBENIUS, &nA));
355       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Approximation error transpose %g (%g / %g, %g)\n", nD / nA, nD, nA, nB));
356       PetscCall(MatDestroy(&D));
357 
358       PetscCall(MatTranspose(Bet, MAT_INPLACE_MATRIX, &Bet));
359       PetscCall(MatAXPY(Be, -1.0, Bet, SAME_NONZERO_PATTERN));
360       PetscCall(MatViewFromOptions(Be, NULL, "-diff_expl_view"));
361       PetscCall(MatNorm(Be, NORM_FROBENIUS, &nB));
362       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Approximation error B - (B^T)^T %g\n", nB));
363       PetscCall(MatDestroy(&Be));
364       PetscCall(MatDestroy(&Bet));
365     }
366   }
367   PetscCall(MatDestroy(&A));
368   PetscCall(MatDestroy(&B));
369   PetscCall(MatDestroy(&C));
370   PetscCall(PetscFree(Cdata));
371   PetscCall(PetscFree(Adata));
372   PetscCall(PetscFinalize());
373   return 0;
374 }
375 
376 /*TEST
377 
378    build:
379      requires: h2opus
380 
381 #tests from kernel
382    test:
383      requires: h2opus
384      nsize: 1
385      suffix: 1
386      args: -n {{17 33}} -kernel 1 -dim {{1 2 3}} -symm {{0 1}} -checkexpl -bgpu 0
387 
388    test:
389      requires: h2opus
390      nsize: 1
391      suffix: 1_ld
392      output_file: output/ex66_1.out
393      args: -n 33 -kernel 1 -dim 1 -lda 13 -ldc 11 -ldu 17 -symm 0 -checkexpl -bgpu 0
394 
395    test:
396      requires: h2opus cuda
397      nsize: 1
398      suffix: 1_cuda
399      output_file: output/ex66_1.out
400      args: -n {{17 33}} -kernel 1 -dim {{1 2 3}} -symm {{0 1}} -checkexpl -bgpu 1
401 
402    test:
403      requires: h2opus cuda
404      nsize: 1
405      suffix: 1_cuda_ld
406      output_file: output/ex66_1.out
407      args: -n 33 -kernel 1 -dim 1 -lda 13 -ldc 11 -ldu 17 -symm 0 -checkexpl -bgpu 1
408 
409    test:
410      requires: h2opus
411      nsize: {{2 3}}
412      suffix: 1_par
413      args: -n 64 -symm -kernel 1 -dim 1 -ldc 12 -testlayout {{0 1}} -bgpu 0 -cgpu 0
414 
415    test:
416      requires: h2opus cuda
417      nsize: {{2 3}}
418      suffix: 1_par_cuda
419      args: -n 64 -symm -kernel 1 -dim 1 -ldc 12 -testlayout {{0 1}} -bgpu {{0 1}} -cgpu {{0 1}}
420      output_file: output/ex66_1_par.out
421 
422 #tests from matrix sampling (parallel or unsymmetric not supported)
423    test:
424      requires: h2opus
425      nsize: 1
426      suffix: 2
427      args: -n {{17 33}} -kernel 0 -dim 2 -symm 1 -checkexpl -bgpu 0
428 
429    test:
430      requires: h2opus cuda
431      nsize: 1
432      suffix: 2_cuda
433      output_file: output/ex66_2.out
434      args: -n {{17 33}} -kernel 0 -dim 2 -symm 1 -checkexpl -bgpu {{0 1}} -agpu {{0 1}}
435 
436 #tests view operation
437    test:
438      requires: h2opus !cuda
439      filter: grep -v " MPI process" | grep -v "\[" | grep -v "\]"
440      nsize: {{1 2 3}}
441      suffix: view
442      args: -ng 64 -kernel 1 -dim 2 -symm 1 -checkexpl -B_view -mat_h2opus_leafsize 17 -mat_h2opus_normsamples 13 -mat_h2opus_indexmap_view ::ascii_matlab -mat_approximate_norm_samples 2 -mat_h2opus_normsamples 2
443 
444    test:
445      requires: h2opus cuda
446      filter: grep -v " MPI process" | grep -v "\[" | grep -v "\]"
447      nsize: {{1 2 3}}
448      suffix: view_cuda
449      args: -ng 64 -kernel 1 -dim 2 -symm 1 -checkexpl -bgpu -B_view -mat_h2opus_leafsize 17 -mat_h2opus_normsamples 13 -mat_h2opus_indexmap_view ::ascii_matlab -mat_approximate_norm_samples 2 -mat_h2opus_normsamples 2
450 
451 TEST*/
452