cannam@127: /* cannam@127: * Copyright (c) 2003, 2007-14 Matteo Frigo cannam@127: * Copyright (c) 2003, 2007-14 Massachusetts Institute of Technology cannam@127: * cannam@127: * This program is free software; you can redistribute it and/or modify cannam@127: * it under the terms of the GNU General Public License as published by cannam@127: * the Free Software Foundation; either version 2 of the License, or cannam@127: * (at your option) any later version. cannam@127: * cannam@127: * This program is distributed in the hope that it will be useful, cannam@127: * but WITHOUT ANY WARRANTY; without even the implied warranty of cannam@127: * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the cannam@127: * GNU General Public License for more details. cannam@127: * cannam@127: * You should have received a copy of the GNU General Public License cannam@127: * along with this program; if not, write to the Free Software cannam@127: * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA cannam@127: * cannam@127: */ cannam@127: cannam@127: cannam@127: #include "verify.h" cannam@127: cannam@127: /* copy real A into real B, using output stride of A and input stride of B */ cannam@127: typedef struct { cannam@127: dotens2_closure k; cannam@127: R *ra; cannam@127: R *rb; cannam@127: } cpyr_closure; cannam@127: cannam@127: static void cpyr0(dotens2_closure *k_, cannam@127: int indxa, int ondxa, int indxb, int ondxb) cannam@127: { cannam@127: cpyr_closure *k = (cpyr_closure *)k_; cannam@127: k->rb[indxb] = k->ra[ondxa]; cannam@127: UNUSED(indxa); UNUSED(ondxb); cannam@127: } cannam@127: cannam@127: static void cpyr(R *ra, const bench_tensor *sza, cannam@127: R *rb, const bench_tensor *szb) cannam@127: { cannam@127: cpyr_closure k; cannam@127: k.k.apply = cpyr0; cannam@127: k.ra = ra; k.rb = rb; cannam@127: bench_dotens2(sza, szb, &k.k); cannam@127: } cannam@127: cannam@127: /* copy unpacked halfcomplex A[n] into packed-complex B[n], using output stride cannam@127: of A and input stride of B. Only copies non-redundant half; other cannam@127: half must be copied via mkhermitian. */ cannam@127: typedef struct { cannam@127: dotens2_closure k; cannam@127: int n; cannam@127: int as; cannam@127: int scalea; cannam@127: R *ra, *ia; cannam@127: R *rb, *ib; cannam@127: } cpyhc2_closure; cannam@127: cannam@127: static void cpyhc20(dotens2_closure *k_, cannam@127: int indxa, int ondxa, int indxb, int ondxb) cannam@127: { cannam@127: cpyhc2_closure *k = (cpyhc2_closure *)k_; cannam@127: int i, n = k->n; cannam@127: int scalea = k->scalea; cannam@127: int as = k->as * scalea; cannam@127: R *ra = k->ra + ondxa * scalea, *ia = k->ia + ondxa * scalea; cannam@127: R *rb = k->rb + indxb, *ib = k->ib + indxb; cannam@127: UNUSED(indxa); UNUSED(ondxb); cannam@127: cannam@127: for (i = 0; i < n/2 + 1; ++i) { cannam@127: rb[2*i] = ra[as*i]; cannam@127: ib[2*i] = ia[as*i]; cannam@127: } cannam@127: } cannam@127: cannam@127: static void cpyhc2(R *ra, R *ia, cannam@127: const bench_tensor *sza, const bench_tensor *vecsza, cannam@127: int scalea, cannam@127: R *rb, R *ib, const bench_tensor *szb) cannam@127: { cannam@127: cpyhc2_closure k; cannam@127: BENCH_ASSERT(sza->rnk <= 1); cannam@127: k.k.apply = cpyhc20; cannam@127: k.n = tensor_sz(sza); cannam@127: k.scalea = scalea; cannam@127: if (!BENCH_FINITE_RNK(sza->rnk) || sza->rnk == 0) cannam@127: k.as = 0; cannam@127: else cannam@127: k.as = sza->dims[0].os; cannam@127: k.ra = ra; k.ia = ia; k.rb = rb; k.ib = ib; cannam@127: bench_dotens2(vecsza, szb, &k.k); cannam@127: } cannam@127: cannam@127: /* icpyhc2 is the inverse of cpyhc2 */ cannam@127: cannam@127: static void icpyhc20(dotens2_closure *k_, cannam@127: int indxa, int ondxa, int indxb, int ondxb) cannam@127: { cannam@127: cpyhc2_closure *k = (cpyhc2_closure *)k_; cannam@127: int i, n = k->n; cannam@127: int scalea = k->scalea; cannam@127: int as = k->as * scalea; cannam@127: R *ra = k->ra + indxa * scalea, *ia = k->ia + indxa * scalea; cannam@127: R *rb = k->rb + ondxb, *ib = k->ib + ondxb; cannam@127: UNUSED(ondxa); UNUSED(indxb); cannam@127: cannam@127: for (i = 0; i < n/2 + 1; ++i) { cannam@127: ra[as*i] = rb[2*i]; cannam@127: ia[as*i] = ib[2*i]; cannam@127: } cannam@127: } cannam@127: cannam@127: static void icpyhc2(R *ra, R *ia, cannam@127: const bench_tensor *sza, const bench_tensor *vecsza, cannam@127: int scalea, cannam@127: R *rb, R *ib, const bench_tensor *szb) cannam@127: { cannam@127: cpyhc2_closure k; cannam@127: BENCH_ASSERT(sza->rnk <= 1); cannam@127: k.k.apply = icpyhc20; cannam@127: k.n = tensor_sz(sza); cannam@127: k.scalea = scalea; cannam@127: if (!BENCH_FINITE_RNK(sza->rnk) || sza->rnk == 0) cannam@127: k.as = 0; cannam@127: else cannam@127: k.as = sza->dims[0].is; cannam@127: k.ra = ra; k.ia = ia; k.rb = rb; k.ib = ib; cannam@127: bench_dotens2(vecsza, szb, &k.k); cannam@127: } cannam@127: cannam@127: typedef struct { cannam@127: dofft_closure k; cannam@127: bench_problem *p; cannam@127: } dofft_rdft2_closure; cannam@127: cannam@127: static void rdft2_apply(dofft_closure *k_, cannam@127: bench_complex *in, bench_complex *out) cannam@127: { cannam@127: dofft_rdft2_closure *k = (dofft_rdft2_closure *)k_; cannam@127: bench_problem *p = k->p; cannam@127: bench_tensor *totalsz, *pckdsz, *totalsz_swap, *pckdsz_swap; cannam@127: bench_tensor *probsz2, *totalsz2, *pckdsz2; cannam@127: bench_tensor *probsz2_swap, *totalsz2_swap, *pckdsz2_swap; cannam@127: bench_real *ri, *ii, *ro, *io; cannam@127: int n2, totalscale; cannam@127: cannam@127: totalsz = tensor_append(p->vecsz, p->sz); cannam@127: pckdsz = verify_pack(totalsz, 2); cannam@127: n2 = tensor_sz(totalsz); cannam@127: if (BENCH_FINITE_RNK(p->sz->rnk) && p->sz->rnk > 0) cannam@127: n2 = (n2 / p->sz->dims[p->sz->rnk - 1].n) * cannam@127: (p->sz->dims[p->sz->rnk - 1].n / 2 + 1); cannam@127: ri = (bench_real *) p->in; cannam@127: ro = (bench_real *) p->out; cannam@127: cannam@127: if (BENCH_FINITE_RNK(p->sz->rnk) && p->sz->rnk > 0 && n2 > 0) { cannam@127: probsz2 = tensor_copy_sub(p->sz, p->sz->rnk - 1, 1); cannam@127: totalsz2 = tensor_copy_sub(totalsz, 0, totalsz->rnk - 1); cannam@127: pckdsz2 = tensor_copy_sub(pckdsz, 0, pckdsz->rnk - 1); cannam@127: } cannam@127: else { cannam@127: probsz2 = mktensor(0); cannam@127: totalsz2 = tensor_copy(totalsz); cannam@127: pckdsz2 = tensor_copy(pckdsz); cannam@127: } cannam@127: cannam@127: totalsz_swap = tensor_copy_swapio(totalsz); cannam@127: pckdsz_swap = tensor_copy_swapio(pckdsz); cannam@127: totalsz2_swap = tensor_copy_swapio(totalsz2); cannam@127: pckdsz2_swap = tensor_copy_swapio(pckdsz2); cannam@127: probsz2_swap = tensor_copy_swapio(probsz2); cannam@127: cannam@127: /* confusion: the stride is the distance between complex elements cannam@127: when using interleaved format, but it is the distance between cannam@127: real elements when using split format */ cannam@127: if (p->split) { cannam@127: ii = p->ini ? (bench_real *) p->ini : ri + n2; cannam@127: io = p->outi ? (bench_real *) p->outi : ro + n2; cannam@127: totalscale = 1; cannam@127: } else { cannam@127: ii = p->ini ? (bench_real *) p->ini : ri + 1; cannam@127: io = p->outi ? (bench_real *) p->outi : ro + 1; cannam@127: totalscale = 2; cannam@127: } cannam@127: cannam@127: if (p->sign < 0) { /* R2HC */ cannam@127: int N, vN, i; cannam@127: cpyr(&c_re(in[0]), pckdsz, ri, totalsz); cannam@127: after_problem_rcopy_from(p, ri); cannam@127: doit(1, p); cannam@127: after_problem_hccopy_to(p, ro, io); cannam@127: if (k->k.recopy_input) cannam@127: cpyr(ri, totalsz_swap, &c_re(in[0]), pckdsz_swap); cannam@127: cpyhc2(ro, io, probsz2, totalsz2, totalscale, cannam@127: &c_re(out[0]), &c_im(out[0]), pckdsz2); cannam@127: N = tensor_sz(p->sz); cannam@127: vN = tensor_sz(p->vecsz); cannam@127: for (i = 0; i < vN; ++i) cannam@127: mkhermitian(out + i*N, p->sz->rnk, p->sz->dims, 1); cannam@127: } cannam@127: else { /* HC2R */ cannam@127: icpyhc2(ri, ii, probsz2, totalsz2, totalscale, cannam@127: &c_re(in[0]), &c_im(in[0]), pckdsz2); cannam@127: after_problem_hccopy_from(p, ri, ii); cannam@127: doit(1, p); cannam@127: after_problem_rcopy_to(p, ro); cannam@127: if (k->k.recopy_input) cannam@127: cpyhc2(ri, ii, probsz2_swap, totalsz2_swap, totalscale, cannam@127: &c_re(in[0]), &c_im(in[0]), pckdsz2_swap); cannam@127: mkreal(out, tensor_sz(pckdsz)); cannam@127: cpyr(ro, totalsz, &c_re(out[0]), pckdsz); cannam@127: } cannam@127: cannam@127: tensor_destroy(totalsz); cannam@127: tensor_destroy(pckdsz); cannam@127: tensor_destroy(totalsz_swap); cannam@127: tensor_destroy(pckdsz_swap); cannam@127: tensor_destroy(probsz2); cannam@127: tensor_destroy(totalsz2); cannam@127: tensor_destroy(pckdsz2); cannam@127: tensor_destroy(probsz2_swap); cannam@127: tensor_destroy(totalsz2_swap); cannam@127: tensor_destroy(pckdsz2_swap); cannam@127: } cannam@127: cannam@127: void verify_rdft2(bench_problem *p, int rounds, double tol, errors *e) cannam@127: { cannam@127: C *inA, *inB, *inC, *outA, *outB, *outC, *tmp; cannam@127: int n, vecn, N; cannam@127: dofft_rdft2_closure k; cannam@127: cannam@127: BENCH_ASSERT(p->kind == PROBLEM_REAL); cannam@127: cannam@127: if (!BENCH_FINITE_RNK(p->sz->rnk) || !BENCH_FINITE_RNK(p->vecsz->rnk)) cannam@127: return; /* give up */ cannam@127: cannam@127: k.k.apply = rdft2_apply; cannam@127: k.k.recopy_input = 0; cannam@127: k.p = p; cannam@127: cannam@127: if (rounds == 0) cannam@127: rounds = 20; /* default value */ cannam@127: cannam@127: n = tensor_sz(p->sz); cannam@127: vecn = tensor_sz(p->vecsz); cannam@127: N = n * vecn; cannam@127: cannam@127: inA = (C *) bench_malloc(N * sizeof(C)); cannam@127: inB = (C *) bench_malloc(N * sizeof(C)); cannam@127: inC = (C *) bench_malloc(N * sizeof(C)); cannam@127: outA = (C *) bench_malloc(N * sizeof(C)); cannam@127: outB = (C *) bench_malloc(N * sizeof(C)); cannam@127: outC = (C *) bench_malloc(N * sizeof(C)); cannam@127: tmp = (C *) bench_malloc(N * sizeof(C)); cannam@127: cannam@127: e->i = impulse(&k.k, n, vecn, inA, inB, inC, outA, outB, outC, cannam@127: tmp, rounds, tol); cannam@127: e->l = linear(&k.k, 1, N, inA, inB, inC, outA, outB, outC, cannam@127: tmp, rounds, tol); cannam@127: cannam@127: e->s = 0.0; cannam@127: if (p->sign < 0) cannam@127: e->s = dmax(e->s, tf_shift(&k.k, 1, p->sz, n, vecn, p->sign, cannam@127: inA, inB, outA, outB, cannam@127: tmp, rounds, tol, TIME_SHIFT)); cannam@127: else cannam@127: e->s = dmax(e->s, tf_shift(&k.k, 1, p->sz, n, vecn, p->sign, cannam@127: inA, inB, outA, outB, cannam@127: tmp, rounds, tol, FREQ_SHIFT)); cannam@127: cannam@127: if (!p->in_place && !p->destroy_input) cannam@127: preserves_input(&k.k, p->sign < 0 ? mkreal : mkhermitian1, cannam@127: N, inA, inB, outB, rounds); cannam@127: cannam@127: bench_free(tmp); cannam@127: bench_free(outC); cannam@127: bench_free(outB); cannam@127: bench_free(outA); cannam@127: bench_free(inC); cannam@127: bench_free(inB); cannam@127: bench_free(inA); cannam@127: } cannam@127: cannam@127: void accuracy_rdft2(bench_problem *p, int rounds, int impulse_rounds, cannam@127: double t[6]) cannam@127: { cannam@127: dofft_rdft2_closure k; cannam@127: int n; cannam@127: C *a, *b; cannam@127: cannam@127: BENCH_ASSERT(p->kind == PROBLEM_REAL); cannam@127: BENCH_ASSERT(p->sz->rnk == 1); cannam@127: BENCH_ASSERT(p->vecsz->rnk == 0); cannam@127: cannam@127: k.k.apply = rdft2_apply; cannam@127: k.k.recopy_input = 0; cannam@127: k.p = p; cannam@127: n = tensor_sz(p->sz); cannam@127: cannam@127: a = (C *) bench_malloc(n * sizeof(C)); cannam@127: b = (C *) bench_malloc(n * sizeof(C)); cannam@127: accuracy_test(&k.k, p->sign < 0 ? mkreal : mkhermitian1, p->sign, cannam@127: n, a, b, rounds, impulse_rounds, t); cannam@127: bench_free(b); cannam@127: bench_free(a); cannam@127: }