Homec4science
Diffusion tamaas (features/materials) Browse

tamaas/src/types.cpp
87f2d0f8cfb1features/materials

Last Change
Raw File

types.cpp
View Options

/**
*
* @author Lucas Frérot <lucas.frerot@epfl.ch>
*
* @section LICENSE
*
* Copyright (©) 2016 EPFL (Ecole Polytechnique Fédérale de
* Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des
* Solides)
*
* Tamaas 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.
*
* Tamaas 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 Lesser General Public License for more
* details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Tamaas. If not, see <http://www.gnu.org/licenses/>.
*
*/
/* -------------------------------------------------------------------------- */
#include "types.hh"
#include <cmath>
/* -------------------------------------------------------------------------- */
__BEGIN_TAMAAS__
/* -------------------------------------------------------------------------- */
/* VectorProxy definitions */
/* -------------------------------------------------------------------------- */
template <typename T>
VectorProxy<T>::VectorProxy(T * mem, UInt n) : Proxy<T, 1>() {
this->data.wrapMemory(mem, n);
this->n[0] = n;
}
/* -------------------------------------------------------------------------- */
template <typename T>
VectorProxy<T>::~VectorProxy() {}
/* -------------------------------------------------------------------------- */
template <typename T>
const T & VectorProxy<T>::operator()(UInt i) const {
return Grid<T, 1>::operator()(i, 0);
}
/* -------------------------------------------------------------------------- */
template <typename T>
T & VectorProxy<T>::operator()(UInt i) {
return Grid<T, 1>::operator()(i, 0);
}
/* -------------------------------------------------------------------------- */
template <typename T>
void VectorProxy<T>::mul(const MatrixProxy<T> &mat,
const VectorProxy<T> &vec) {
const UInt * n_mat = mat.sizes();
const UInt * n_vec = vec.sizes();
TAMAAS_ASSERT(n_mat[1] == n_vec[0], "matrix and vector sizes do not match");
TAMAAS_ASSERT(n_mat[0] == this->n[0],
"matrix and result vector sizes do not match");
// No parallelism here
for (UInt i = 0 ; i < n_mat[0] ; ++i) {
for (UInt j = 0 ; j < n_vec[0] ; ++j) {
(*this)(i) += mat(i, j)*vec(j);
}
}
}
/* -------------------------------------------------------------------------- */
template <typename T>
T VectorProxy<T>::dot(const VectorProxy<T> &vec) {
TAMAAS_ASSERT(this->n[0] == vec.sizes()[0], "vector sizes do not match");
T res = 0;
// No parallelism here
for (UInt i = 0 ; i < this->n[0] ; i++)
res += (*this)(i) * vec(i);
return res;
}
/* -------------------------------------------------------------------------- */
template <typename T>
T VectorProxy<T>::l2norm() const {
T res = 0;
// No parallelism here
for (UInt i = 0 ; i < this->n[0] ; i++)
res += (*this)(i)*(*this)(i);
return std::sqrt(res);
}
/* -------------------------------------------------------------------------- */
/* MatrixProxy definitions */
/* -------------------------------------------------------------------------- */
template <typename T>
MatrixProxy<T>::MatrixProxy(T * mem, UInt n, UInt m) : Proxy<T, 2>() {
this->data.wrapMemory(mem, n*m);
this->n[0] = n;
this->n[1] = m;
}
/* -------------------------------------------------------------------------- */
template <typename T>
MatrixProxy<T>::~MatrixProxy() {}
/* -------------------------------------------------------------------------- */
template <typename T>
const T & MatrixProxy<T>::operator()(UInt i, UInt j) const {
return Grid<T, 2>::operator()(i, j, 0);
}
/* -------------------------------------------------------------------------- */
template <typename T>
T & MatrixProxy<T>::operator()(UInt i, UInt j) {
return Grid<T, 2>::operator()(i, j, 0);
}
/* -------------------------------------------------------------------------- */
#define TYPE_INSTANCIATE(type) template class VectorProxy<type>; \
template class MatrixProxy<type>
TYPE_INSTANCIATE(Real);
TYPE_INSTANCIATE(UInt);
TYPE_INSTANCIATE(Complex);
TYPE_INSTANCIATE(int);
TYPE_INSTANCIATE(unsigned long);
#undef TYPE_INSTANCIATE
__END_TAMAAS__
c4science · Help