Skip to content
Snippets Groups Projects
Commit 57bd4a70 authored by patavirt's avatar patavirt
Browse files

Adjust to const + rename SimpleVector

parent 14b4732f
No related branches found
No related tags found
No related merge requests found
...@@ -15,15 +15,15 @@ ...@@ -15,15 +15,15 @@
#include "array.hpp" #include "array.hpp"
template <typename Scalar> template <typename Scalar, typename ScalarW = typename std::remove_const<Scalar>::type>
inline Matrix<Scalar,4,4> inline Matrix<ScalarW,4,4>
get_Q_matrix(array::ArrayView<Scalar, array::Shape<2,2,2> > Q) get_Q_matrix(array::ArrayView<Scalar, array::Shape<2,2,2> > Q)
{ {
auto g = Q.part(0u).matrix(); auto g = Q.part(0u).matrix();
auto gt = Q.part(1u).matrix(); auto gt = Q.part(1u).matrix();
Matrix<double,2,2> I; Matrix<double,2,2> I;
Matrix<Scalar,2,2> N, Nt; Matrix<ScalarW,2,2> N, Nt;
Matrix<Scalar,4,4> Qm; Matrix<ScalarW,4,4> Qm;
I << 1, 0, I << 1, 0,
0, 1; 0, 1;
...@@ -39,13 +39,13 @@ get_Q_matrix(array::ArrayView<Scalar, array::Shape<2,2,2> > Q) ...@@ -39,13 +39,13 @@ get_Q_matrix(array::ArrayView<Scalar, array::Shape<2,2,2> > Q)
return Qm; return Qm;
} }
template <typename Scalar, typename Shape, typename Shape2> template <typename Scalar, typename Shape, typename Shape2, typename ScalarW = typename std::remove_const<Scalar>::type>
inline Scalar S_2(array::ArrayView<Scalar, Shape> Q, array::ArrayView<Complex, Shape2> U, double h) inline ScalarW S_2(array::ArrayView<Scalar, Shape> Q, array::ArrayView<const Complex, Shape2> U, double h)
{ {
const int N_i[4] = {0, 1, -1, 0}; const int N_i[4] = {0, 1, -1, 0};
const int N_j[4] = {1, 0, 0, -1}; const int N_j[4] = {1, 0, 0, -1};
const size_t invk[4] = {3, 2, 1, 0}; const size_t invk[4] = {3, 2, 1, 0};
Scalar S = Scalar(0); ScalarW S = Scalar(0);
size_t nx = Q.dim(0); size_t nx = Q.dim(0);
size_t ny = Q.dim(1); size_t ny = Q.dim(1);
......
...@@ -16,7 +16,7 @@ ...@@ -16,7 +16,7 @@
#include <cppad/example/cppad_eigen.hpp> #include <cppad/example/cppad_eigen.hpp>
#include "array.hpp" #include "array.hpp"
#include "simplevector.hpp" #include "vector.hpp"
typedef std::complex<double> Complex; typedef std::complex<double> Complex;
typedef CppAD::AD<double> ADDouble; typedef CppAD::AD<double> ADDouble;
......
...@@ -14,8 +14,8 @@ using namespace array::python; ...@@ -14,8 +14,8 @@ using namespace array::python;
class CLASS_Action class CLASS_Action
{ {
public: public:
typedef PyArrayView<Complex,Shape<Dynamic,Dynamic,4,4,4> > UType; typedef PyArrayView<const Complex,Shape<Dynamic,Dynamic,4,4,4> > UType;
typedef PyArrayView<Complex,Shape<Dynamic,Dynamic,2,2,2> > QType; typedef PyArrayView<const Complex,Shape<Dynamic,Dynamic,2,2,2> > QType;
private: private:
UType U_; UType U_;
...@@ -39,9 +39,9 @@ public: ...@@ -39,9 +39,9 @@ public:
if (Qval.dim(0) != nx_ || Qval.dim(1) != ny_) if (Qval.dim(0) != nx_ || Qval.dim(1) != ny_)
throw std::out_of_range("Q array wrong size"); throw std::out_of_range("Q array wrong size");
StoredComplexArray<ADComplex, Shape<Dynamic,Dynamic,2,2,2> > Array<ADComplex, Shape<Dynamic,Dynamic,2,2,2>, std::vector<ADDouble> >
Q({nx_, ny_, 2, 2, 2}); Q({nx_, ny_, 2, 2, 2});
StoredComplexArray<ADComplex, Shape<1> > res; Array<ADComplex, Shape<1>, std::vector<ADDouble> > res;
auto a = Qval.reshape<Dynamic>({Qval.size()}); auto a = Qval.reshape<Dynamic>({Qval.size()});
auto b = Q.reshape<Dynamic>({Q.size()}); auto b = Q.reshape<Dynamic>({Q.size()});
...@@ -57,14 +57,13 @@ public: ...@@ -57,14 +57,13 @@ public:
auto grad(QType Q) auto grad(QType Q)
{ {
SimpleVector<double> x(reinterpret_cast<double *>(Q.data()), Vector<double> x(reinterpret_cast<double *>(const_cast<Complex *>(Q.data())), 2 * Q.size());
2 * Q.size());
CppAD::ADFun<double> f = eval_ad(Q); CppAD::ADFun<double> f = eval_ad(Q);
auto jac = f.Jacobian(x); auto jac = f.Jacobian(x);
std::array<size_t, 1> shape = {jac.size()/2}; std::array<size_t, 1> shape = {jac.size()/2};
return StoredArray<Complex, Shape<Dynamic>, 0, true, SimpleVector<double> >(std::move(jac), shape); return Array<Complex, Shape<Dynamic>, Vector<double> >(std::move(jac), shape);
} }
}; };
......
...@@ -25,10 +25,10 @@ void dump(ArrayView<Scalar,Shape> Q) ...@@ -25,10 +25,10 @@ void dump(ArrayView<Scalar,Shape> Q)
int main0() int main0()
{ {
StoredComplexArray<Complex,Shape<2,2,2,2,2> > Qval; Array<Complex,Shape<2,2,2,2,2>,std::vector<double> > Qval;
StoredComplexArray<ADComplex,Shape<2,2,2,2,2> > Q; Array<ADComplex,Shape<2,2,2,2,2>,std::vector<ADDouble> > Q;
StoredArray<Complex,Shape<2,2,4,4,4> > U; Array<Complex,Shape<2,2,4,4,4> > U;
for (size_t i = 0; i < Qval.dim(0); ++i) { for (size_t i = 0; i < Qval.dim(0); ++i) {
for (size_t j = 0; i < Qval.dim(1); ++i) { for (size_t j = 0; i < Qval.dim(1); ++i) {
...@@ -51,8 +51,10 @@ int main0() ...@@ -51,8 +51,10 @@ int main0()
CppAD::Independent(Q.storage()); CppAD::Independent(Q.storage());
StoredComplexArray<ADComplex,Shape<1> > res; Array<ADComplex,Shape<1>,std::vector<ADDouble> > res;
res(0u) = S_2(Q, U, 0.1);
ArrayView<const Complex,Shape<2,2,4,4,4> > U_const = U;
res(0u) = S_2(Q, U_const, 0.1);
CppAD::ADFun<double> f(Q.storage(), res.storage()); CppAD::ADFun<double> f(Q.storage(), res.storage());
...@@ -70,7 +72,7 @@ int main0() ...@@ -70,7 +72,7 @@ int main0()
int main() int main()
{ {
StoredArray<ADDouble,Shape<8,8> > Q; Array<ADDouble,Shape<8,8> > Q;
for (size_t i = 0; i < Q.dim(0); ++i) { for (size_t i = 0; i < Q.dim(0); ++i) {
for (size_t j = 0; i < Q.dim(1); ++i) { for (size_t j = 0; i < Q.dim(1); ++i) {
...@@ -80,7 +82,7 @@ int main() ...@@ -80,7 +82,7 @@ int main()
CppAD::Independent(Q.storage()); CppAD::Independent(Q.storage());
StoredArray<ADDouble,Shape<1> > res; Array<ADDouble,Shape<1> > res;
res(0u) = laplacian(Q); res(0u) = laplacian(Q);
CppAD::ADFun<double> f(Q.storage(), res.storage()); CppAD::ADFun<double> f(Q.storage(), res.storage());
......
...@@ -10,8 +10,10 @@ ...@@ -10,8 +10,10 @@
#include <vector> #include <vector>
namespace array {
template <class Scalar> template <class Scalar>
class SimpleVector class Vector
{ {
public: public:
typedef Scalar value_type; typedef Scalar value_type;
...@@ -40,15 +42,15 @@ private: ...@@ -40,15 +42,15 @@ private:
} }
public: public:
SimpleVector() Vector()
: storage_(), data_(storage_.data()), size_(storage_.size()) : storage_(), data_(storage_.data()), size_(storage_.size())
{} {}
SimpleVector(size_t n) Vector(size_t n)
: storage_(n), data_(storage_.data()), size_(storage_.size()) : storage_(n), data_(storage_.data()), size_(storage_.size())
{} {}
SimpleVector(const SimpleVector& other) Vector(const Vector& other)
: storage_(other.storage_), data_(other.data_), size_(other.size_) : storage_(other.storage_), data_(other.data_), size_(other.size_)
{ {
if (other.own_data()) { if (other.own_data()) {
...@@ -59,18 +61,18 @@ public: ...@@ -59,18 +61,18 @@ public:
} }
} }
SimpleVector(SimpleVector&& other) Vector(Vector&& other)
: storage_(std::move(other.storage_)), data_(other.data_), size_(other.size_) : storage_(std::move(other.storage_)), data_(other.data_), size_(other.size_)
{ {
other.data_ = other.storage_.data(); other.data_ = other.storage_.data();
other.size_ = other.storage_.size(); other.size_ = other.storage_.size();
} }
SimpleVector(Scalar *data, size_t size) Vector(Scalar *data, size_t size)
: storage_(), data_(data), size_(size) : storage_(), data_(data), size_(size)
{} {}
SimpleVector& operator=(const SimpleVector& other) Vector& operator=(const Vector& other)
{ {
if (other.own_data()) { if (other.own_data()) {
storage_ = other.storage_; storage_ = other.storage_;
...@@ -84,7 +86,7 @@ public: ...@@ -84,7 +86,7 @@ public:
return *this; return *this;
} }
SimpleVector& operator=(SimpleVector&& other) Vector& operator=(Vector&& other)
{ {
if (other.own_data()) { if (other.own_data()) {
storage_ = std::move(other.storage_); storage_ = std::move(other.storage_);
...@@ -131,4 +133,6 @@ public: ...@@ -131,4 +133,6 @@ public:
} }
}; };
}
#endif // SIMPLEVECTOR_HPP_ #endif // SIMPLEVECTOR_HPP_
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment