Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Another charge against pointer to pointer #1312

Merged
merged 11 commits into from
Jul 18, 2021
Merged

Another charge against pointer to pointer #1312

Show file tree
Hide file tree
Changes from 7 commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
4f8fccc
create a type for "matrix views" (wraps a pointer) to avoid **
pcarruscag Jun 27, 2021
6971ab0
small fix
pcarruscag Jun 27, 2021
31e94a0
const correctness and less duplication in CNumerics
pcarruscag Jun 27, 2021
7d1eb8c
correct brief
pcarruscag Jun 27, 2021
72e119a
more cleaning + update regressions
pcarruscag Jun 27, 2021
973b331
Merge branch 'gradient_type' of https://github.com/su2code/SU2 into g…
pcarruscag Jun 27, 2021
a488554
not needed member variable
pcarruscag Jun 27, 2021
936fe67
Merge branch 'develop' into gradient_type
pcarruscag Jun 29, 2021
e348065
Merge branch 'develop' into gradient_type
pcarruscag Jul 15, 2021
520b978
Apply suggestions from code review
pcarruscag Jul 15, 2021
0ad7a02
Update SU2_CFD/include/numerics/CNumerics.hpp
pcarruscag Jul 15, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
87 changes: 54 additions & 33 deletions Common/include/containers/container_decorators.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,43 @@

#include "C2DContainer.hpp"

/*!
* \brief Class to represent a matrix (without owning the data, this just wraps a pointer).
*/
template<class T>
class CMatrixView {
public:
using Scalar = typename std::remove_const<T>::type;
using Index = unsigned long;

private:
T* m_ptr;
Index m_cols;

public:
CMatrixView(T* ptr = nullptr, Index cols = 0) : m_ptr(ptr), m_cols(cols) {}

template<class U> friend class CMatrixView;
template<class U>
CMatrixView(const CMatrixView<U>& other) : m_ptr(other.m_ptr), m_cols(other.m_cols) {}

explicit CMatrixView(su2matrix<Scalar>& mat) : m_ptr(mat.data()), m_cols(mat.cols()) {}

template<class U = T, su2enable_if<std::is_const<U>::value> = 0>
explicit CMatrixView(const su2matrix<Scalar>& mat) : m_ptr(mat.data()), m_cols(mat.cols()) {}

const Scalar* operator[] (Index i) const noexcept { return &m_ptr[i*m_cols]; }
const Scalar& operator() (Index i, Index j) const noexcept { return m_ptr[i*m_cols + j]; }

template<class U = T, su2enable_if<!std::is_const<U>::value> = 0>
Scalar* operator[] (Index i) noexcept { return &m_ptr[i*m_cols]; }

template<class U = T, su2enable_if<!std::is_const<U>::value> = 0>
Scalar& operator() (Index i, Index j) noexcept { return m_ptr[i*m_cols + j]; }
Comment on lines +58 to +65
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess this detects whether a var on the receiving side is const and then uses the latter versions of [] and () operators?
Edit: Or the matrix view - object has to be created with this property already and then it decides and what to do? I guess the latter

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hmmm no..
If the matrix view object is used in a const context then it returns const refs to su2double. Otherwise it will "try" to return non const ref.
However, if the view was returned by a const container, it cannot do that. Hence the use of sfinae to disable the non const version.


friend CMatrixView operator+ (CMatrixView mv, Index incr) { return CMatrixView(mv[incr], mv.m_cols); }
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So this allows for the +1 logic which is equivalent to (i,j) logic to allow for similar use as with pointers. This of course might break less code of people and their code out there when pulling this, so 👍 ... but on the other hand might give some false security that the underlying datatype is still the same?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Well the type is not really important in this case, just the interface. And there is too much of this +stuff, so less work for me xD.

};

/*!
* \class C3DContainerDecorator
* \brief Decorate a matrix type (Storage) with 3 dimensions.
Expand All @@ -44,6 +81,9 @@ class C3DContainerDecorator {
static constexpr bool IsRowMajor = true;
static constexpr bool IsColumnMajor = false;

using Matrix = CMatrixView<Scalar>;
using ConstMatrix = CMatrixView<const Scalar>;

using CInnerIter = typename Storage::CInnerIter;
template<class T, size_t N>
using CInnerIterGather = typename Storage::template CInnerIterGather<simd::Array<T,N> >;
Expand Down Expand Up @@ -78,6 +118,18 @@ class C3DContainerDecorator {
Scalar& operator() (Index i, Index j, Index k) noexcept { return m_storage(i, j*m_innerSz + k); }
const Scalar& operator() (Index i, Index j, Index k) const noexcept { return m_storage(i, j*m_innerSz + k); }

/*!
* \brief Matrix access.
*/
Matrix operator[] (Index i) noexcept { return Matrix(m_storage[i], m_innerSz); }
ConstMatrix operator[] (Index i) const noexcept { return ConstMatrix(m_storage[i], m_innerSz); }
Comment on lines +124 to +125
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So C3DContainerDecorator owns the data and can return MatrixView objects that provide access to the 2D structure for one entry after the first layer (or the 3rd dimension to stay in this dimension language)?

Copy link
Member Author

@pcarruscag pcarruscag Jul 15, 2021 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think I was influenced by Eigen (inception style) for the "view"


/*!
* \brief Matrix access with an offset.
*/
Matrix operator() (Index i, Index j) noexcept { return Matrix(m_storage[i]+j*m_innerSz, m_innerSz); }
ConstMatrix operator() (Index i, Index j) const noexcept { return ConstMatrix(m_storage[i]+j*m_innerSz, m_innerSz); }

/*!
* \brief Get a scalar iterator to the inner-most dimension of the container.
*/
Expand Down Expand Up @@ -105,42 +157,11 @@ class C3DContainerDecorator {
};

/*!
* \brief Some typedefs for the
* \brief Some typedefs for 3D containers
*/
using C3DIntMatrix = C3DContainerDecorator<su2matrix<unsigned long> >;
using C3DDoubleMatrix = C3DContainerDecorator<su2activematrix>;

/*!
* \class CVectorOfMatrix
* \brief This contrived container is used to store small matrices in a contiguous manner
* but still present the "su2double**" interface to the outside world.
* The "interface" part should be replaced by something more efficient, e.g. a "matrix view".
*/
class CVectorOfMatrix: public C3DDoubleMatrix {
private:
su2matrix<Scalar*> interface;
Comment on lines -115 to -121
Copy link
Member Author

@pcarruscag pcarruscag Jun 27, 2021 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The main motivation for this change is to store "vectors of matrices" e.g. gradients in a contiguous way without needing an extra matrix of pointers to provide the su2double** compatibility with other areas of the code.


public:
CVectorOfMatrix() = default;

CVectorOfMatrix(Index length, Index rows, Index cols, Scalar value = 0) noexcept {
resize(length, rows, cols, value);
}

void resize(Index length, Index rows, Index cols, Scalar value = 0) noexcept {
C3DDoubleMatrix::resize(length, rows, cols, value);
interface.resize(length,rows);
for(Index i=0; i<length; ++i)
for(Index j=0; j<rows; ++j)
interface(i,j) = &(*this)(i,j,0);
}

/*!
* \brief Matrix-wise access.
*/
Scalar** operator[] (Index i) noexcept { return interface[i]; }
const Scalar* const* operator[] (Index i) const noexcept { return interface[i]; }
};
using CVectorOfMatrix = C3DDoubleMatrix;

/*!
* \class C2DDummyLastView
Expand Down
2 changes: 1 addition & 1 deletion Common/include/geometry/dual_grid/CPoint.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -809,7 +809,7 @@ class CPoint {
* \param[in] iPoint - Index of the point.
* \return Grid velocity gradient at the point.
*/
inline const su2double* const* GetGridVel_Grad(unsigned long iPoint) const { return GridVel_Grad[iPoint]; }
inline CMatrixView<const su2double> GetGridVel_Grad(unsigned long iPoint) const { return GridVel_Grad[iPoint]; }

/*!
* \brief Set the adjoint values of the (geometric) coordinates.
Expand Down
Loading