1# Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2# the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3# reserved. See files LICENSE and NOTICE for details. 4# 5# This file is part of CEED, a collection of benchmarks, miniapps, software 6# libraries and APIs for efficient high-order finite element and spectral 7# element discretizations for exascale applications. For more information and 8# source code availability see http://github.com/ceed. 9# 10# The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11# a collaborative effort of two U.S. Department of Energy organizations (Office 12# of Science and the National Nuclear Security Administration) responsible for 13# the planning and preparation of a capable exascale ecosystem, including 14# software, applications, hardware, advanced system engineering and early 15# testbed platforms, in support of the nation's exascale computing imperative. 16 17# @file 18# Test Ceed ElemRestriction functionality 19 20import os 21import libceed 22import numpy as np 23import check 24 25# ------------------------------------------------------------------------------- 26# Test creation, use, and destruction of an element restriction 27# ------------------------------------------------------------------------------- 28 29 30def test_200(ceed_resource): 31 ceed = libceed.Ceed(ceed_resource) 32 33 num_elem = 3 34 35 x = ceed.Vector(num_elem + 1) 36 a = np.arange(10, 10 + num_elem + 1, dtype=ceed.scalar_type()) 37 x.set_array(a, cmode=libceed.USE_POINTER) 38 39 ind = np.zeros(2 * num_elem, dtype="int32") 40 for i in range(num_elem): 41 ind[2 * i + 0] = i 42 ind[2 * i + 1] = i + 1 43 r = ceed.ElemRestriction(num_elem, 2, 1, 1, num_elem + 1, ind, 44 cmode=libceed.USE_POINTER) 45 46 y = ceed.Vector(2 * num_elem) 47 y.set_value(0) 48 49 r.apply(x, y) 50 51 with y.array_read() as y_array: 52 for i in range(2 * num_elem): 53 assert 10 + (i + 1) // 2 == y_array[i] 54 55# ------------------------------------------------------------------------------- 56# Test creation, use, and destruction of a strided element restriction 57# ------------------------------------------------------------------------------- 58 59 60def test_201(ceed_resource): 61 ceed = libceed.Ceed(ceed_resource) 62 63 num_elem = 3 64 65 x = ceed.Vector(2 * num_elem) 66 a = np.arange(10, 10 + 2 * num_elem, dtype=ceed.scalar_type()) 67 x.set_array(a, cmode=libceed.USE_POINTER) 68 69 strides = np.array([1, 2, 2], dtype="int32") 70 r = ceed.StridedElemRestriction(num_elem, 2, 1, 2 * num_elem, strides) 71 72 y = ceed.Vector(2 * num_elem) 73 y.set_value(0) 74 75 r.apply(x, y) 76 77 with y.array_read() as y_array: 78 for i in range(2 * num_elem): 79 assert 10 + i == y_array[i] 80 81# ------------------------------------------------------------------------------- 82# Test creation and destruction of a blocked element restriction 83# ------------------------------------------------------------------------------- 84 85 86def test_202(ceed_resource, capsys): 87 ceed = libceed.Ceed(ceed_resource) 88 89 num_elem = 8 90 elem_size = 2 91 num_blk = 2 92 blk_size = 5 93 94 x = ceed.Vector(num_elem + 1) 95 a = np.arange(10, 10 + num_elem + 1, dtype=ceed.scalar_type()) 96 x.set_array(a, cmode=libceed.COPY_VALUES) 97 98 ind = np.zeros(2 * num_elem, dtype="int32") 99 for i in range(num_elem): 100 ind[2 * i + 0] = i 101 ind[2 * i + 1] = i + 1 102 r = ceed.BlockedElemRestriction(num_elem, elem_size, blk_size, 1, 1, 103 num_elem + 1, ind, cmode=libceed.COPY_VALUES) 104 105 y = ceed.Vector(num_blk * blk_size * elem_size) 106 y.set_value(0) 107 108 # NoTranspose 109 r.apply(x, y) 110 layout = r.get_layout() 111 with y.array_read() as y_array: 112 for i in range(elem_size): 113 for j in range(1): 114 for k in range(num_elem): 115 block = int(k / blk_size) 116 elem = k % blk_size 117 indx = (i * blk_size + elem) * \ 118 layout[0] + j * blk_size * layout[1] + \ 119 block * blk_size * layout[2] 120 assert y_array[indx] == a[ind[k * elem_size + i]] 121 122 x.set_value(0) 123 r.T.apply(y, x) 124 with x.array_read() as x_array: 125 for i in range(num_elem + 1): 126 assert x_array[i] == (10 + i) * (2.0 if i > 127 0 and i < num_elem else 1.0) 128 129# ------------------------------------------------------------------------------- 130# Test creation, use, and destruction of a blocked element restriction 131# ------------------------------------------------------------------------------- 132 133 134def test_208(ceed_resource): 135 ceed = libceed.Ceed(ceed_resource) 136 137 num_elem = 8 138 elem_size = 2 139 num_blk = 2 140 blk_size = 5 141 142 x = ceed.Vector(num_elem + 1) 143 a = np.arange(10, 10 + num_elem + 1, dtype=ceed.scalar_type()) 144 x.set_array(a, cmode=libceed.COPY_VALUES) 145 146 ind = np.zeros(2 * num_elem, dtype="int32") 147 for i in range(num_elem): 148 ind[2 * i + 0] = i 149 ind[2 * i + 1] = i + 1 150 r = ceed.BlockedElemRestriction(num_elem, elem_size, blk_size, 1, 1, 151 num_elem + 1, ind, cmode=libceed.COPY_VALUES) 152 153 y = ceed.Vector(blk_size * elem_size) 154 y.set_value(0) 155 156 # NoTranspose 157 r.apply_block(1, x, y) 158 layout = r.get_layout() 159 with y.array_read() as y_array: 160 for i in range(elem_size): 161 for j in range(1): 162 for k in range(blk_size, num_elem): 163 block = int(k / blk_size) 164 elem = k % blk_size 165 indx = (i * blk_size + elem) * layout[0] + j * blk_size * \ 166 layout[1] + block * blk_size * \ 167 layout[2] - blk_size * elem_size 168 assert y_array[indx] == a[ind[k * elem_size + i]] 169 170 x.set_value(0) 171 r.T.apply_block(1, y, x) 172 with x.array_read() as x_array: 173 for i in range(blk_size, num_elem + 1): 174 assert x_array[i] == (10 + i) * (2.0 if i > 175 blk_size and i < num_elem else 1.0) 176 177# ------------------------------------------------------------------------------- 178# Test getting the multiplicity of the indices in an element restriction 179# ------------------------------------------------------------------------------- 180 181 182def test_209(ceed_resource): 183 ceed = libceed.Ceed(ceed_resource) 184 185 num_elem = 3 186 187 ind = np.zeros(4 * num_elem, dtype="int32") 188 for i in range(num_elem): 189 ind[4 * i + 0] = i * 3 + 0 190 ind[4 * i + 1] = i * 3 + 1 191 ind[4 * i + 2] = i * 3 + 2 192 ind[4 * i + 3] = i * 3 + 3 193 r = ceed.ElemRestriction(num_elem, 4, 1, 1, 3 * num_elem + 1, ind, 194 cmode=libceed.USE_POINTER) 195 196 mult = r.get_multiplicity() 197 198 with mult.array_read() as mult_array: 199 for i in range(3 * num_elem + 1): 200 val = 1 + (1 if (i > 0 and i < 3 * num_elem and i % 3 == 0) else 0) 201 assert val == mult_array[i] 202 203# ------------------------------------------------------------------------------- 204# Test creation and view of an element restriction 205# ------------------------------------------------------------------------------- 206 207 208def test_210(ceed_resource, capsys): 209 ceed = libceed.Ceed(ceed_resource) 210 211 num_elem = 3 212 213 ind = np.zeros(2 * num_elem, dtype="int32") 214 for i in range(num_elem): 215 ind[2 * i + 0] = i + 0 216 ind[2 * i + 1] = i + 1 217 r = ceed.ElemRestriction(num_elem, 2, 1, 1, num_elem + 1, ind, 218 cmode=libceed.USE_POINTER) 219 220 print(r) 221 222 stdout, stderr, ref_stdout = check.output(capsys) 223 assert not stderr 224 assert stdout == ref_stdout 225 226# ------------------------------------------------------------------------------- 227# Test creation and view of a strided element restriction 228# ------------------------------------------------------------------------------- 229 230 231def test_211(ceed_resource, capsys): 232 ceed = libceed.Ceed(ceed_resource) 233 234 num_elem = 3 235 236 strides = np.array([1, 2, 2], dtype="int32") 237 r = ceed.StridedElemRestriction(num_elem, 2, 1, num_elem + 1, strides) 238 239 print(r) 240 241 stdout, stderr, ref_stdout = check.output(capsys) 242 assert not stderr 243 assert stdout == ref_stdout 244 245# ------------------------------------------------------------------------------- 246# Test creation and view of a blocked strided element restriction 247# ------------------------------------------------------------------------------- 248 249 250def test_212(ceed_resource, capsys): 251 ceed = libceed.Ceed(ceed_resource) 252 253 num_elem = 3 254 255 strides = np.array([1, 2, 2], dtype="int32") 256 r = ceed.BlockedStridedElemRestriction( 257 num_elem, 2, 2, 1, num_elem + 1, strides) 258 259 print(r) 260 261 stdout, stderr, ref_stdout = check.output(capsys) 262 assert not stderr 263 assert stdout == ref_stdout 264 265# ------------------------------------------------------------------------------- 266