Basix 0.4.0

Home     Installation     Demos     C++ docs     Python docs

precompute.h
1 // Copyright (c) 2020 Matthew Scroggs
2 // FEniCS Project
3 // SPDX-License-Identifier: MIT
4 
5 #pragma once
6 
7 #include <tuple>
8 #include <vector>
9 #include <xtensor/xtensor.hpp>
10 #include <xtl/xspan.hpp>
11 
14 {
15 
64 std::vector<std::size_t>
65 prepare_permutation(const std::vector<std::size_t>& perm);
66 
118 template <typename E>
119 void apply_permutation(const std::vector<std::size_t>& perm,
120  const xtl::span<E>& data, std::size_t offset = 0,
121  std::size_t block_size = 1)
122 {
123  for (std::size_t b = 0; b < block_size; ++b)
124  {
125  for (std::size_t i = 0; i < perm.size(); ++i)
126  {
127  std::swap(data[block_size * (offset + i) + b],
128  data[block_size * (offset + perm[i]) + b]);
129  }
130  }
131 }
132 
139 template <typename E>
140 void apply_permutation_to_transpose(const std::vector<std::size_t>& perm,
141  const xtl::span<E>& data,
142  std::size_t offset = 0,
143  std::size_t block_size = 1)
144 {
145  const std::size_t dim = perm.size();
146  const std::size_t data_size
147  = (data.size() + (dim < block_size ? block_size - dim : 0)) / block_size;
148  for (std::size_t b = 0; b < block_size; ++b)
149  {
150  for (std::size_t i = 0; i < dim; ++i)
151  {
152  std::swap(data[data_size * b + offset + i],
153  data[data_size * b + offset + perm[i]]);
154  }
155  }
156 }
157 
251 std::tuple<std::vector<std::size_t>, std::vector<double>,
252  xt::xtensor<double, 2>>
253 prepare_matrix(const xt::xtensor<double, 2>& matrix);
254 
318 template <typename T, typename E>
319 void apply_matrix(const std::tuple<std::vector<std::size_t>, std::vector<T>,
320  xt::xtensor<T, 2>>& matrix,
321  const xtl::span<E>& data, std::size_t offset = 0,
322  std::size_t block_size = 1)
323 {
324  const std::vector<std::size_t>& v_size_t = std::get<0>(matrix);
325  const std::vector<T>& v_t = std::get<1>(matrix);
326  const xt::xtensor<T, 2>& M = std::get<2>(matrix);
327 
328  const std::size_t dim = v_size_t.size();
329  apply_permutation(v_size_t, data, offset, block_size);
330  for (std::size_t b = 0; b < block_size; ++b)
331  {
332  for (std::size_t i = 0; i < dim; ++i)
333  {
334  data[block_size * (offset + i) + b] *= v_t[i];
335  for (std::size_t j = 0; j < dim; ++j)
336  {
337  data[block_size * (offset + i) + b]
338  += M(i, j) * data[block_size * (offset + j) + b];
339  }
340  }
341  }
342 }
343 
350 template <typename T, typename E>
352  const std::tuple<std::vector<std::size_t>, std::vector<T>,
353  xt::xtensor<T, 2>>& matrix,
354  const xtl::span<E>& data, std::size_t offset = 0,
355  std::size_t block_size = 1)
356 {
357  const std::vector<std::size_t>& v_size_t = std::get<0>(matrix);
358  const std::vector<T>& v_t = std::get<1>(matrix);
359  const xt::xtensor<T, 2>& M = std::get<2>(matrix);
360 
361  const std::size_t dim = v_size_t.size();
362  const std::size_t data_size
363  = (data.size() + (dim < block_size ? block_size - dim : 0)) / block_size;
364  apply_permutation_to_transpose(v_size_t, data, offset, block_size);
365  for (std::size_t b = 0; b < block_size; ++b)
366  {
367  for (std::size_t i = 0; i < dim; ++i)
368  {
369  data[data_size * b + offset + i] *= v_t[i];
370  for (std::size_t j = 0; j < dim; ++j)
371  {
372  data[data_size * b + offset + i]
373  += M(i, j) * data[data_size * b + offset + j];
374  }
375  }
376  }
377 }
378 
379 } // namespace basix::precompute
basix::precompute::prepare_permutation
std::vector< std::size_t > prepare_permutation(const std::vector< std::size_t > &perm)
Definition: precompute.cpp:13
basix::precompute::apply_permutation_to_transpose
void apply_permutation_to_transpose(const std::vector< std::size_t > &perm, const xtl::span< E > &data, std::size_t offset=0, std::size_t block_size=1)
Definition: precompute.h:140
basix::precompute::prepare_matrix
std::tuple< std::vector< std::size_t >, std::vector< double >, xt::xtensor< double, 2 > > prepare_matrix(const xt::xtensor< double, 2 > &matrix)
Definition: precompute.cpp:28
basix::precompute::apply_permutation
void apply_permutation(const std::vector< std::size_t > &perm, const xtl::span< E > &data, std::size_t offset=0, std::size_t block_size=1)
Definition: precompute.h:119
basix::precompute::apply_matrix_to_transpose
void apply_matrix_to_transpose(const std::tuple< std::vector< std::size_t >, std::vector< T >, xt::xtensor< T, 2 >> &matrix, const xtl::span< E > &data, std::size_t offset=0, std::size_t block_size=1)
Definition: precompute.h:351
basix::precompute::apply_matrix
void apply_matrix(const std::tuple< std::vector< std::size_t >, std::vector< T >, xt::xtensor< T, 2 >> &matrix, const xtl::span< E > &data, std::size_t offset=0, std::size_t block_size=1)
Definition: precompute.h:319
basix::precompute
Matrix and permutation precomputation.
Definition: precompute.h:13