using int instead of unsigned int in all template parameters

fdaPDE/core.h

@@ -14,8 +14,8 @@
 // You should have received a copy of the GNU General Public License
 // along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
-#ifndef __FDAPDE_CORE_MOUDLE_H__
-#define __FDAPDE_CORE_MOUDLE_H__
+#ifndef __FDAPDE_CORE_MODULE_H__
+#define __FDAPDE_CORE_MODULE_H__
 
 #include "utils.h"
 #include "fields.h"
@@ -26,4 +26,4 @@
 #include "splines.h"
 #include "multithreading.h"
 
-#endif   // __FDAPDE_CORE_MOUDLE_H__
+#endif   // __FDAPDE_CORE_MODULE_H__

fdaPDE/fields.h

@@ -14,8 +14,8 @@
 // You should have received a copy of the GNU General Public License
 // along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
-#ifndef __FDAPDE_FIELDS_MOUDLE_H__
-#define __FDAPDE_FIELDS_MOUDLE_H__
+#ifndef __FDAPDE_FIELDS_MODULE_H__
+#define __FDAPDE_FIELDS_MODULE_H__
 
 #include "fields/differentiable_field.h"
 #include "fields/dot_product.h"
@@ -27,4 +27,4 @@
 #include "fields/vector_expressions.h"
 #include "fields/vector_field.h"
 
-#endif   // __FDAPDE_FIELDS_MOUDLE_H__
+#endif   // __FDAPDE_FIELDS_MODULE_H__

fdaPDE/finite_elements/basis/lagrangian_element.h

@@ -27,17 +27,17 @@ namespace fdapde {
 namespace core {
 
 // A Lagrangian basis of degree R over an M-dimensional element
-template <unsigned int M_, unsigned int R_> class LagrangianElement {
+template <int M_, int R_> class LagrangianElement {
    private:
     std::array<std::array<double, M_>, ct_binomial_coefficient(M_ + R_, R_)> nodes_;   // nodes of the Lagrangian basis
     std::array<MultivariatePolynomial<M_, R_>, ct_binomial_coefficient(M_ + R_, R_)> basis_;
     // solves the Vandermonde system for the computation of polynomial coefficients
     void compute_coefficients_(const std::array<std::array<double, M_>, ct_binomial_coefficient(M_ + R_, R_)>& nodes);
    public:
     // compile time informations
-    static constexpr unsigned int R = R_; // basis order
-    static constexpr unsigned int M = M_; // input space dimension
-    static constexpr unsigned int n_basis = ct_binomial_coefficient(M + R, R);
+    static constexpr int R = R_; // basis order
+    static constexpr int M = M_; // input space dimension
+    static constexpr int n_basis = ct_binomial_coefficient(M + R, R);
     typedef MultivariatePolynomial<M, R> ElementType;
     typedef typename std::array<MultivariatePolynomial<M, R>, n_basis>::const_iterator const_iterator;
 
@@ -59,11 +59,11 @@ template <unsigned int M_, unsigned int R_> class LagrangianElement {
 // implementative details
 
 // compute coefficients via Vandermonde matrix
-template <unsigned int M_, unsigned int R_>
+template <int M_, int R_>
 void LagrangianElement<M_, R_>::compute_coefficients_(
   const std::array<std::array<double, M_>, ct_binomial_coefficient(M_ + R_, R_)>& nodes) {
     // build vandermonde matrix
-    constexpr unsigned int n_basis = ct_binomial_coefficient(M_ + R_, R_);
+    constexpr int n_basis = ct_binomial_coefficient(M_ + R_, R_);
     constexpr std::array<std::array<unsigned, M_>, n_basis> poly_table = MultivariatePolynomial<M, R>::poly_table;
 
     // Vandermonde matrix construction

fdaPDE/finite_elements/basis/multivariate_polynomial.h

@@ -47,10 +47,9 @@ namespace core {
 //   2   0   0  | 2
 //
 // ct_poly_exp() computes the above table for any choice of N and R at compile time.
-template <unsigned int N, unsigned int R>
-using PolyTable = std::array<std::array<unsigned, N>, ct_binomial_coefficient(R + N, R)>;
+template <int N, int R> using PolyTable = std::array<std::array<unsigned, N>, ct_binomial_coefficient(R + N, R)>;
 
-template <unsigned int N, unsigned int R> constexpr PolyTable<N, R> ct_poly_exp() {
+template <int N, int R> constexpr PolyTable<N, R> ct_poly_exp() {
     const int monomials = ct_binomial_coefficient(R + N, R);   // number of monomials in polynomial p(x)
     // initialize empty array
     PolyTable<N, R> coeff {};
@@ -98,10 +97,10 @@ template <unsigned int N, unsigned int R> constexpr PolyTable<N, R> ct_poly_exp(
 //     PolyTable                         PolyGradTable
 //
 // ct_grad_exp() computes the PolyGradTable above from a PolyTable for any choice of N and R at compile time.
-template <unsigned int N, unsigned int R>
+template <int N, int R>
 using PolyGradTable = std::array<std::array<std::array<unsigned, N>, ct_binomial_coefficient(R + N, R)>, N>;
 
-template <unsigned int N, unsigned int R> constexpr PolyGradTable<N, R> ct_grad_exp(const PolyTable<N, R> poly_table) {
+template <int N, int R> constexpr PolyGradTable<N, R> ct_grad_exp(const PolyTable<N, R> poly_table) {
     const int monomials = ct_binomial_coefficient(R + N, R);   // number of monomials in polynomial p(x)
     // initialize empty array
     PolyGradTable<N,R> coeff {};
@@ -120,9 +119,9 @@ template <unsigned int N, unsigned int R> constexpr PolyGradTable<N, R> ct_grad_
 }
 
 // recursive template based unfolding of monomial product x1^i1*x2^i2*...*xN^iN.
-template <unsigned int N,   // template recursion loop variable
-	  typename P,       // point where to evaluate the monomial
-	  typename V>       // a row of the polynomial PolyTable (i.e. an array of coefficients [i1 i2 ... iN])
+template <int N,        // template recursion loop variable
+          typename P,   // point where to evaluate the monomial
+          typename V>   // a row of the polynomial PolyTable (i.e. an array of coefficients [i1 i2 ... iN])
 struct MonomialProduct {
     static constexpr double unfold(const P& p, const V& v) {
         return v[N] == 0 ? MonomialProduct<N - 1, P, V>::unfold(p, v) :
@@ -135,26 +134,26 @@ template <typename P, typename V> struct MonomialProduct<0, P, V> {   // base ca
 };
 
 // unfold the sum of monomials at compile time to produce the complete polynomial expression
-template <unsigned int I,   // template recursion loop variable
-	  unsigned int N,   // total number of monomials to unfold
-	  unsigned int M,   // polynomial space dimension
-	  typename P,       // point where to evaluate the polynomial
-	  typename V>       // the whole polynomial PolyTable
+template <int I,        // template recursion loop variable
+          int N,        // total number of monomials to unfold
+          int M,        // polynomial space dimension
+          typename P,   // point where to evaluate the polynomial
+          typename V>   // the whole polynomial PolyTable
 struct MonomialSum {
     static constexpr double unfold(const std::array<double, N>& c, const P& p, const V& v) {
         return (c[I] * MonomialProduct<M - 1, SVector<M>, std::array<unsigned, M>>::unfold(p, v[I])) +
                MonomialSum<I - 1, N, M, P, V>::unfold(c, p, v);
     }
 };
 // end of recursion
-template <unsigned int N, unsigned int M, typename P, typename V> struct MonomialSum<0, N, M, P, V> {
+template <int N, int M, typename P, typename V> struct MonomialSum<0, N, M, P, V> {
     static constexpr double unfold(const std::array<double, N>& c, const P& p, const V& v) {
         return c[0] * MonomialProduct<M - 1, SVector<M>, std::array<unsigned, M>>::unfold(p, v[0]);
     }
 };
 
 // functor implementing the derivative of a multivariate N-dimensional polynomial of degree R along a given direction
-template <unsigned int N, unsigned int R>
+template <int N, int R>
 class PolynomialDerivative : public VectorExpr<N, N, PolynomialDerivative<N, R>> {
    private:
     // compile time informations
@@ -184,25 +183,25 @@ class PolynomialDerivative : public VectorExpr<N, N, PolynomialDerivative<N, R>>
 };
 
 // class representing a multivariate polynomial of degree R defined over a space of dimension N
-template <unsigned int N, unsigned int R>
+template <int N, int R>
 class MultivariatePolynomial : public ScalarExpr<N, MultivariatePolynomial<N, R>> {
    private:
-    static const constexpr unsigned MON = ct_binomial_coefficient(R + N, R);
+    static const constexpr int MON = ct_binomial_coefficient(R + N, R);
     std::array<double, MON> coeff_vector_;   // vector of coefficients
     VectorField<N, N, PolynomialDerivative<N, R>> gradient_;
    public:
     // compute this at compile time once, let public access
     static constexpr PolyTable<N, R> poly_table = ct_poly_exp<N, R>();
-    static constexpr unsigned int order = R;                   // polynomial order
-    static constexpr unsigned int input_space_dimension = N;   // input space dimension
+    static constexpr int order = R;                   // polynomial order
+    static constexpr int input_space_dimension = N;   // input space dimension
 
     // constructor
     MultivariatePolynomial() = default;
     MultivariatePolynomial(const std::array<double, MON>& coeff_vector) : coeff_vector_(coeff_vector) {
         // prepare gradient vector
         std::array<PolynomialDerivative<N, R>, N> gradient;
         // define i-th element of gradient field
-        for (size_t i = 0; i < N; ++i) { gradient[i] = PolynomialDerivative<N, R>(coeff_vector, i); }
+        for (std::size_t i = 0; i < N; ++i) { gradient[i] = PolynomialDerivative<N, R>(coeff_vector, i); }
         gradient_ = VectorField<N, N, PolynomialDerivative<N, R>>(gradient);
     };
     // evaluate polynomial at point

fdaPDE/linear_algebra.h

@@ -14,12 +14,12 @@
 // You should have received a copy of the GNU General Public License
 // along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
-#ifndef __FDAPDE_LINEAR_ALGEBRA_MOUDLE_H__
-#define __FDAPDE_LINEAR_ALGEBRA_MOUDLE_H__
+#ifndef __FDAPDE_LINEAR_ALGEBRA_MODULE_H__
+#define __FDAPDE_LINEAR_ALGEBRA_MODULE_H__
 
 #include "linear_algebra/kronecker_product.h"
 #include "linear_algebra/smw.h"
 #include "linear_algebra/sparse_block_matrix.h"
 #include "linear_algebra/vector_space.h"
 
-#endif   // __FDAPDE_LINEAR_ALGEBRA_MOUDLE_H__
+#endif   // __FDAPDE_LINEAR_ALGEBRA_MODULE_H__

fdaPDE/linear_algebra/vector_space.h

@@ -27,7 +27,7 @@ namespace core {
 
 // a template class to perform geometric operations in general vector and affine spaces.
 // M is the vector space dimension, N is the dimension of the embedding space
-template <unsigned int M, unsigned int N> class VectorSpace {
+template <int M, int N> class VectorSpace {
    private:
     std::array<SVector<N>, M> basis_ {};   // the set of vectors generating the vector space
     SVector<N> offset_ {};   // a point throught which the vector space passes. (offset_ \neq 0 \implies affine space)
@@ -51,12 +51,12 @@ template <unsigned int M, unsigned int N> class VectorSpace {
 
 // implementation details
 
-template <unsigned int M, unsigned int N>
+template <int M, int N>
 SVector<N> VectorSpace<M, N>::orthogonal_project(const SVector<N>& u, const SVector<N>& v) const {
     return ((v.dot(u) / u.squaredNorm()) * (u.array())).matrix();
 }
 
-template <unsigned int M, unsigned int N> void VectorSpace<M, N>::orthonormalize() {
+template <int M, int N> void VectorSpace<M, N>::orthonormalize() {
     std::array<SVector<N>, M> orthonormal_basis;
     orthonormal_basis[0] = basis_[0] / basis_[0].norm();
     // implementation of the modified Gram-Schmidt method, see theory for details
@@ -70,7 +70,7 @@ template <unsigned int M, unsigned int N> void VectorSpace<M, N>::orthonormalize
 }
 
 // projects the point x onto the space.
-template <unsigned int M, unsigned int N> SVector<M> VectorSpace<M, N>::project_onto(const SVector<N>& x) {
+template <int M, int N> SVector<M> VectorSpace<M, N>::project_onto(const SVector<N>& x) {
     if constexpr (M == N) {   // you are projecting over the same space!
         return x;
     } else {
@@ -82,7 +82,7 @@ template <unsigned int M, unsigned int N> SVector<M> VectorSpace<M, N>::project_
 }
 
 // project the point x into the space
-template <unsigned int M, unsigned int N> SVector<N> VectorSpace<M, N>::project_into(const SVector<N>& x) {
+template <int M, int N> SVector<N> VectorSpace<M, N>::project_into(const SVector<N>& x) {
     // build the projection operator on the space spanned by the basis
     Eigen::Matrix<double, N, Eigen::Dynamic> A;
     A.resize(N, basis_.size());
@@ -93,12 +93,12 @@ template <unsigned int M, unsigned int N> SVector<N> VectorSpace<M, N>::project_
 }
 
 // euclidean distance from x
-template <unsigned int M, unsigned int N> double VectorSpace<M, N>::distance(const SVector<N>& x) {
+template <int M, int N> double VectorSpace<M, N>::distance(const SVector<N>& x) {
     return (x - project_into(x)).squaredNorm();
 }
 
 // develops the linear combination of basis vectors with respect to the given vector of coefficients.
-template <unsigned int M, unsigned int N>
+template <int M, int N>
 SVector<N> VectorSpace<M, N>::operator()(const std::array<double, M>& coeffs) const {
     SVector<N> result = offset_;
     for (std::size_t i = 0; i < M; ++i) result += coeffs[i] * basis_[i];

fdaPDE/mesh.h

@@ -14,8 +14,8 @@
 // You should have received a copy of the GNU General Public License
 // along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
-#ifndef __FDAPDE_MESH_MOUDLE_H__
-#define __FDAPDE_MESH_MOUDLE_H__
+#ifndef __FDAPDE_MESH_MODULE_H__
+#define __FDAPDE_MESH_MODULE_H__
 
 #include "mesh/mesh.h"
 #include "mesh/element.h"
@@ -25,4 +25,4 @@
 #include "mesh/point_location/barycentric_walk.h"
 #include "mesh/point_location/adt.h"
 
-#endif   // __FDAPDE_MESH_MOUDLE_H__
+#endif   // __FDAPDE_MESH_MODULE_H__

fdaPDE/mesh/element.h

@@ -29,24 +29,24 @@ namespace fdapde {
 namespace core {
 
 // M local dimension, N embedding dimension, R order of the element (defaulted to linear finite elements)
-template <unsigned int M, unsigned int N, unsigned int R = 1> class Element;   // forward declaration
+template <int M, int N, int R = 1> class Element;   // forward declaration
 
 // using declarations for manifold specialization
-template <unsigned int R> using SurfaceElement = Element<2, 3, R>;
-template <unsigned int R> using NetworkElement = Element<1, 2, R>;
+template <int R> using SurfaceElement = Element<2, 3, R>;
+template <int R> using NetworkElement = Element<1, 2, R>;
 
 // number of degrees of freedom associated to an M-dimensional element of degree R
-constexpr unsigned int ct_nnodes(const unsigned int M, const unsigned int R) {
+constexpr int ct_nnodes(const int M, const int R) {
     return ct_factorial(M + R) / (ct_factorial(M) * ct_factorial(R));
 }
 // number of vertices of an M-dimensional element
-constexpr unsigned int ct_nvertices(const unsigned int M) { return M + 1; }
+constexpr int ct_nvertices(const int M) { return M + 1; }
 // number of edges of an M-dimensional element
-constexpr unsigned int ct_nedges(const unsigned int M) { return (M * (M + 1)) / 2; }
+constexpr int ct_nedges(const int M) { return (M * (M + 1)) / 2; }
 
 // A single mesh element. This object represents the main **geometrical** abstraction of a physical element.
 // No functional information is carried by instances of this object.
-template <unsigned int M, unsigned int N, unsigned int R> class Element {
+template <int M, int N, int R> class Element {
    private:
     std::size_t ID_;                                         // ID of this element
     std::array<std::size_t, ct_nvertices(M)> node_ids_ {};   // ID of nodes composing the element
@@ -60,11 +60,11 @@ template <unsigned int M, unsigned int N, unsigned int R> class Element {
     SMatrix<N, M> J_;       // [J_]_ij = (coords_[j][i] - coords_[0][i])
     SMatrix<M, N> inv_J_;   // J^{-1} (Penrose pseudo-inverse for manifold elements)
    public:
-    static constexpr unsigned int nodes = ct_nnodes(M, R);
-    static constexpr unsigned int vertices = ct_nvertices(M);
-    static constexpr unsigned int local_dimension = M;
-    static constexpr unsigned int embedding_dimension = N;
-    static constexpr unsigned int order = R;
+    static constexpr int nodes = ct_nnodes(M, R);
+    static constexpr int vertices = ct_nvertices(M);
+    static constexpr int local_dimension = M;
+    static constexpr int embedding_dimension = N;
+    static constexpr int order = R;
 
     // constructor
     Element() = default;
@@ -74,7 +74,7 @@ template <unsigned int M, unsigned int N, unsigned int R> class Element {
     // getters
     const std::array<SVector<N>, ct_nvertices(M)> coords() const { return coords_; }
     const std::vector<int> neighbors() const { return neighbors_; }
-    const unsigned int ID() const { return ID_; }
+    const int ID() const { return ID_; }
     Eigen::Matrix<double, N, M> barycentric_matrix() const { return J_; }
     Eigen::Matrix<double, M, N> inv_barycentric_matrix() const { return inv_J_; }
     std::array<std::size_t, ct_nvertices(M)> node_ids() const { return node_ids_; }
@@ -99,7 +99,7 @@ template <unsigned int M, unsigned int N, unsigned int R> class Element {
 
 // implementation details
 
-template <unsigned int M, unsigned int N, unsigned int R>
+template <int M, int N, int R>
 Element<M, N, R>::Element(
   std::size_t ID, const std::array<std::size_t, ct_nvertices(M)>& node_ids,
   const std::array<SVector<N>, ct_nvertices(M)>& coords, const std::vector<int>& neighbors, bool boundary) :
@@ -127,7 +127,7 @@ Element<M, N, R>::Element(
     }
 }
 
-template <unsigned int M, unsigned int N, unsigned int R>
+template <int M, int N, int R>
 SVector<M + 1> Element<M, N, R>::to_barycentric_coords(const SVector<N>& x) const {
     // solve: barycenitrc_matrix_*z = (x-ref)
     SVector<M> z = inv_J_ * (x - coords_[0]);
@@ -137,15 +137,15 @@ SVector<M + 1> Element<M, N, R>::to_barycentric_coords(const SVector<N>& x) cons
     return result;
 }
 
-template <unsigned int M, unsigned int N, unsigned int R> SVector<N> Element<M, N, R>::mid_point() const {
+template <int M, int N, int R> SVector<N> Element<M, N, R>::mid_point() const {
     // The center of gravity of an element has all its barycentric coordinates equal to 1/(M+1).
     // The midpoint of an element can be computed by mapping it to a cartesian reference system
     SVector<M> barycentric_mid_point;
     barycentric_mid_point.fill(1.0 / (M + 1));   // avoid implicit conversion to int
     return J_ * barycentric_mid_point + coords_[0];
 }
 
-template <unsigned int M, unsigned int N, unsigned int R>
+template <int M, int N, int R>
 std::pair<SVector<N>, SVector<N>> Element<M, N, R>::bounding_box() const {
     // define lower-left and upper-right corner of bounding box
     SVector<N> ll, ur;
@@ -164,22 +164,22 @@ std::pair<SVector<N>, SVector<N>> Element<M, N, R>::bounding_box() const {
 }
 
 // check if a point is contained in the element
-template <unsigned int M, unsigned int N, unsigned int R>
+template <int M, int N, int R>
 template <bool is_manifold>
 typename std::enable_if<!is_manifold, bool>::type Element<M, N, R>::contains(const SVector<N>& x) const {
     // A point is inside the element if all its barycentric coordinates are positive
     return (to_barycentric_coords(x).array() >= -10 * std::numeric_limits<double>::epsilon()).all();
 }
 
-template <unsigned int M, unsigned int N, unsigned int R> VectorSpace<M, N> Element<M, N, R>::spanned_space() const {
+template <int M, int N, int R> VectorSpace<M, N> Element<M, N, R>::spanned_space() const {
     // build a basis for the space containing the element
     std::array<SVector<N>, M> basis;
     for (size_t i = 0; i < M; ++i) { basis[i] = (coords_[i + 1] - coords_[0]); }
     return VectorSpace<M, N>(basis, coords_[0]);
 }
 
 // specialization for manifold elements of contains() routine.
-template <unsigned int M, unsigned int N, unsigned int R>
+template <int M, int N, int R>
 template <bool is_manifold>
 typename std::enable_if<is_manifold, bool>::type Element<M, N, R>::contains(const SVector<N>& x) const {
     // check if the point is contained in the affine space spanned by the element