From 0426e8a3755ebb9e1d6994c989959d7e440e0381 Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Mon, 27 Nov 2023 17:49:51 -0500 Subject: [PATCH 01/12] fix assignment for nullptr var_value and for assigning expressions --- stan/math/rev/core/arena_matrix.hpp | 4 ++ stan/math/rev/core/var.hpp | 64 +++++++++++++++++++++++++--- test/unit/math/rev/core/var_test.cpp | 46 ++++++++++++++++++++ 3 files changed, 107 insertions(+), 7 deletions(-) diff --git a/stan/math/rev/core/arena_matrix.hpp b/stan/math/rev/core/arena_matrix.hpp index 3d3b095fc74..b65bed948cd 100644 --- a/stan/math/rev/core/arena_matrix.hpp +++ b/stan/math/rev/core/arena_matrix.hpp @@ -128,6 +128,10 @@ class arena_matrix : public Eigen::Map { Base::operator=(a); return *this; } + template + void hard_copy(const T& x) { + Base::operator=(x); + } }; } // namespace math diff --git a/stan/math/rev/core/var.hpp b/stan/math/rev/core/var.hpp index 4e510bacde0..5cae18abdad 100644 --- a/stan/math/rev/core/var.hpp +++ b/stan/math/rev/core/var.hpp @@ -390,6 +390,8 @@ class var_value> { reverse_pass_callback( [this_vi = this->vi_, other_vi = other.vi_]() mutable { other_vi->adj_ += this_vi->adj_; + // + this_vi->adj_.setZero(); }); } @@ -1020,9 +1022,9 @@ class var_value> { * @param other the value to assign * @return this */ - template * = nullptr, - require_all_plain_type_t* = nullptr, - require_not_same_t, plain_type_t>* = nullptr> + template * = nullptr, + require_all_plain_type_t* = nullptr, + require_not_same_t, plain_type_t>* = nullptr> inline var_value& operator=(const var_value& other) { static_assert( EIGEN_PREDICATE_SAME_MATRIX_SIZE(T, S), @@ -1032,16 +1034,63 @@ class var_value> { } /** - * Assignment of another var value, when either this or the other one does not + * Assignment of another var value, when the `this` does not * contain a plain type. - * @tparam S type of the value in the `var_value` to assing + * @tparam S type of the value in the `var_value` to assign + * @param other the value to assign + * @return this + */ + template * = nullptr, + require_not_plain_type_t* = nullptr, + require_plain_type_t* = nullptr> + inline var_value& operator=(const var_value& other) { + // If vi_ is nullptr then the var needs initialized via copy constructor + if (!(this->vi_)) { + *this = var_value(other); + return *this; + } + arena_t> prev_val(vi_->val_.rows(), vi_->val_.cols()); + prev_val.hard_copy(vi_->val_); + vi_->val_.hard_copy(other.val()); + // no need to change any adjoints - these are just zeros before the reverse + // pass + + reverse_pass_callback( + [this_vi = this->vi_, other_vi = other.vi_, prev_val]() mutable { + this_vi->val_.hard_copy(prev_val); + + // we have no way of detecting aliasing between this->vi_->adj_ and + // other.vi_->adj_, so we must copy adjoint before reseting to zero + + // we can reuse prev_val instead of allocating a new matrix + prev_val.hard_copy(this_vi->adj_); + this_vi->adj_.setZero(); + other_vi->adj_ += prev_val; + }); + return *this; + } + /** + * Assignment of another var value, when either both `this` or other does not + * contain a plain type. + * @tparam S type of the value in the `var_value` to assign * @param other the value to assign * @return this */ template * = nullptr, - require_any_not_plain_type_t* = nullptr> + require_any_not_plain_type_t* = nullptr, + require_not_plain_type_t* = nullptr> inline var_value& operator=(const var_value& other) { + // If vi_ is nullptr then the var needs initialized via copy constructor + if (!(this->vi_)) { + []() STAN_COLD_PATH { + throw std::domain_error( + "var_value::operator=(var_value):" + " Internal Bug! Please report this with an example" + " of your model to the Stan math github repository."); + }(); + } arena_t> prev_val = vi_->val_; vi_->val_ = other.val(); // no need to change any adjoints - these are just zeros before the reverse @@ -1055,13 +1104,14 @@ class var_value> { // other.vi_->adj_, so we must copy adjoint before reseting to zero // we can reuse prev_val instead of allocating a new matrix - prev_val = this_vi->adj_; + prev_val.hard_copy(this_vi->adj_); this_vi->adj_.setZero(); other_vi->adj_ += prev_val; }); return *this; } + /** * No-op to match with Eigen methods which call eval */ diff --git a/test/unit/math/rev/core/var_test.cpp b/test/unit/math/rev/core/var_test.cpp index 3e6576e033c..9bd36e2de73 100644 --- a/test/unit/math/rev/core/var_test.cpp +++ b/test/unit/math/rev/core/var_test.cpp @@ -910,3 +910,49 @@ TEST_F(AgradRev, matrix_compile_time_conversions) { EXPECT_MATRIX_FLOAT_EQ(colvec.val(), rowvec.val()); EXPECT_MATRIX_FLOAT_EQ(x11.val(), rowvec.val()); } + +TEST_F(AgradRev, assign_nan) { + using stan::math::var_value; + using var_vector = var_value>; + using stan::math::var; + Eigen::VectorXd x_val(10); + for (int i = 0; i < 10; ++i) { + x_val(i) = i + 0.1; + } + var_vector x(x_val); + var_vector y = var_vector(Eigen::Matrix::Constant(10, std::numeric_limits::quiet_NaN())); + y = stan::math::head(x, 10); + var sigma = 1.0; + var lp = stan::math::normal_lpdf(y, 0, sigma); + lp.grad(); + Eigen::VectorXd x_ans_adj(10); + for (int i = 0; i < 10; ++i) { + x_ans_adj(i) = -(i + 0.1); + } + EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); + Eigen::VectorXd y_ans_adj = Eigen::VectorXd::Zero(10); + EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); +} + +TEST_F(AgradRev, assign_nullptr_vari) { + using stan::math::var_value; + using var_vector = var_value>; + using stan::math::var; + Eigen::VectorXd x_val(10); + for (int i = 0; i < 10; ++i) { + x_val(i) = i + 0.1; + } + var_vector x(x_val); + var_vector y; + y = stan::math::head(x, 10); + var sigma = 1.0; + var lp = stan::math::normal_lpdf(y, 0, sigma); + lp.grad(); + Eigen::VectorXd x_ans_adj(10); + for (int i = 0; i < 10; ++i) { + x_ans_adj(i) = -(i + 0.1); + } + EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); + Eigen::VectorXd y_ans_adj = Eigen::VectorXd::Zero(10); + EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); +} From f5d49a69932e69077dbbeeecf1c79678dca55e2e Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Mon, 27 Nov 2023 17:59:02 -0500 Subject: [PATCH 02/12] [Jenkins] auto-formatting by clang-format version 10.0.0-4ubuntu1 --- stan/math/rev/core/var.hpp | 14 ++--- test/unit/math/rev/core/var_test.cpp | 81 ++++++++++++++-------------- 2 files changed, 48 insertions(+), 47 deletions(-) diff --git a/stan/math/rev/core/var.hpp b/stan/math/rev/core/var.hpp index 5cae18abdad..f02fbaec05c 100644 --- a/stan/math/rev/core/var.hpp +++ b/stan/math/rev/core/var.hpp @@ -1022,7 +1022,8 @@ class var_value> { * @param other the value to assign * @return this */ - template * = nullptr, + template * = nullptr, require_all_plain_type_t* = nullptr, require_not_same_t, plain_type_t>* = nullptr> inline var_value& operator=(const var_value& other) { @@ -1049,7 +1050,7 @@ class var_value> { if (!(this->vi_)) { *this = var_value(other); return *this; - } + } arena_t> prev_val(vi_->val_.rows(), vi_->val_.cols()); prev_val.hard_copy(vi_->val_); vi_->val_.hard_copy(other.val()); @@ -1086,11 +1087,11 @@ class var_value> { if (!(this->vi_)) { []() STAN_COLD_PATH { throw std::domain_error( - "var_value::operator=(var_value):" - " Internal Bug! Please report this with an example" - " of your model to the Stan math github repository."); + "var_value::operator=(var_value):" + " Internal Bug! Please report this with an example" + " of your model to the Stan math github repository."); }(); - } + } arena_t> prev_val = vi_->val_; vi_->val_ = other.val(); // no need to change any adjoints - these are just zeros before the reverse @@ -1111,7 +1112,6 @@ class var_value> { return *this; } - /** * No-op to match with Eigen methods which call eval */ diff --git a/test/unit/math/rev/core/var_test.cpp b/test/unit/math/rev/core/var_test.cpp index 9bd36e2de73..9cb60ba5e80 100644 --- a/test/unit/math/rev/core/var_test.cpp +++ b/test/unit/math/rev/core/var_test.cpp @@ -912,47 +912,48 @@ TEST_F(AgradRev, matrix_compile_time_conversions) { } TEST_F(AgradRev, assign_nan) { - using stan::math::var_value; - using var_vector = var_value>; - using stan::math::var; - Eigen::VectorXd x_val(10); - for (int i = 0; i < 10; ++i) { - x_val(i) = i + 0.1; - } - var_vector x(x_val); - var_vector y = var_vector(Eigen::Matrix::Constant(10, std::numeric_limits::quiet_NaN())); - y = stan::math::head(x, 10); - var sigma = 1.0; - var lp = stan::math::normal_lpdf(y, 0, sigma); - lp.grad(); - Eigen::VectorXd x_ans_adj(10); - for (int i = 0; i < 10; ++i) { - x_ans_adj(i) = -(i + 0.1); - } - EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); - Eigen::VectorXd y_ans_adj = Eigen::VectorXd::Zero(10); - EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); + using stan::math::var_value; + using var_vector = var_value>; + using stan::math::var; + Eigen::VectorXd x_val(10); + for (int i = 0; i < 10; ++i) { + x_val(i) = i + 0.1; + } + var_vector x(x_val); + var_vector y = var_vector(Eigen::Matrix::Constant( + 10, std::numeric_limits::quiet_NaN())); + y = stan::math::head(x, 10); + var sigma = 1.0; + var lp = stan::math::normal_lpdf(y, 0, sigma); + lp.grad(); + Eigen::VectorXd x_ans_adj(10); + for (int i = 0; i < 10; ++i) { + x_ans_adj(i) = -(i + 0.1); + } + EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); + Eigen::VectorXd y_ans_adj = Eigen::VectorXd::Zero(10); + EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); } TEST_F(AgradRev, assign_nullptr_vari) { - using stan::math::var_value; - using var_vector = var_value>; - using stan::math::var; - Eigen::VectorXd x_val(10); - for (int i = 0; i < 10; ++i) { - x_val(i) = i + 0.1; - } - var_vector x(x_val); - var_vector y; - y = stan::math::head(x, 10); - var sigma = 1.0; - var lp = stan::math::normal_lpdf(y, 0, sigma); - lp.grad(); - Eigen::VectorXd x_ans_adj(10); - for (int i = 0; i < 10; ++i) { - x_ans_adj(i) = -(i + 0.1); - } - EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); - Eigen::VectorXd y_ans_adj = Eigen::VectorXd::Zero(10); - EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); + using stan::math::var_value; + using var_vector = var_value>; + using stan::math::var; + Eigen::VectorXd x_val(10); + for (int i = 0; i < 10; ++i) { + x_val(i) = i + 0.1; + } + var_vector x(x_val); + var_vector y; + y = stan::math::head(x, 10); + var sigma = 1.0; + var lp = stan::math::normal_lpdf(y, 0, sigma); + lp.grad(); + Eigen::VectorXd x_ans_adj(10); + for (int i = 0; i < 10; ++i) { + x_ans_adj(i) = -(i + 0.1); + } + EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); + Eigen::VectorXd y_ans_adj = Eigen::VectorXd::Zero(10); + EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); } From a8ab023e2470e29737d1ad3fd50bb2d154a4ee81 Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Mon, 27 Nov 2023 18:02:07 -0500 Subject: [PATCH 03/12] update docs --- stan/math/rev/core/arena_matrix.hpp | 5 +++++ stan/math/rev/core/var.hpp | 2 +- 2 files changed, 6 insertions(+), 1 deletion(-) diff --git a/stan/math/rev/core/arena_matrix.hpp b/stan/math/rev/core/arena_matrix.hpp index b65bed948cd..a0b24725fb7 100644 --- a/stan/math/rev/core/arena_matrix.hpp +++ b/stan/math/rev/core/arena_matrix.hpp @@ -128,6 +128,11 @@ class arena_matrix : public Eigen::Map { Base::operator=(a); return *this; } + /** + * Forces hard copying matrices into an arena matrix + * @tparam T Any type assignable to `Base` + * @param x the values to write to `this` + */ template void hard_copy(const T& x) { Base::operator=(x); diff --git a/stan/math/rev/core/var.hpp b/stan/math/rev/core/var.hpp index f02fbaec05c..5553296383e 100644 --- a/stan/math/rev/core/var.hpp +++ b/stan/math/rev/core/var.hpp @@ -390,7 +390,7 @@ class var_value> { reverse_pass_callback( [this_vi = this->vi_, other_vi = other.vi_]() mutable { other_vi->adj_ += this_vi->adj_; - // + // Reset the adjoint for `this` to replicate SoA before assignment this_vi->adj_.setZero(); }); } From 0d2a46a0732ba86e4569814a3f74579c6dfd7ce8 Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Mon, 27 Nov 2023 18:03:07 -0500 Subject: [PATCH 04/12] [Jenkins] auto-formatting by clang-format version 10.0.0-4ubuntu1 --- stan/math/rev/core/arena_matrix.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/stan/math/rev/core/arena_matrix.hpp b/stan/math/rev/core/arena_matrix.hpp index a0b24725fb7..d6a83f2b6f6 100644 --- a/stan/math/rev/core/arena_matrix.hpp +++ b/stan/math/rev/core/arena_matrix.hpp @@ -129,7 +129,7 @@ class arena_matrix : public Eigen::Map { return *this; } /** - * Forces hard copying matrices into an arena matrix + * Forces hard copying matrices into an arena matrix * @tparam T Any type assignable to `Base` * @param x the values to write to `this` */ From b1eb578777f3107f1f08cdd4908f8f7eccad00eb Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Tue, 28 Nov 2023 11:42:03 -0500 Subject: [PATCH 05/12] add setZero() to matrix_cl types --- stan/math/opencl/matrix_cl.hpp | 31 +++++++++++++++++++++++++++++++ 1 file changed, 31 insertions(+) diff --git a/stan/math/opencl/matrix_cl.hpp b/stan/math/opencl/matrix_cl.hpp index bf276bb19f5..bbb963634c6 100644 --- a/stan/math/opencl/matrix_cl.hpp +++ b/stan/math/opencl/matrix_cl.hpp @@ -501,6 +501,37 @@ class matrix_cl : public matrix_cl_base { */ ~matrix_cl() { wait_for_read_write_events(); } + /** + * Set the values of a `matrix_cl` to zero. + */ + void setZero() { + if (this->size() == 0) { + return; + } + cl::Event zero_event; + this->wait_for_read_write_events(); + const std::size_t write_events_size = this->write_events().size(); + const std::size_t read_events_size = this->read_events().size(); + const std::size_t read_write_size = write_events_size + read_events_size; + std::vector read_write_events(read_write_size, cl::Event{}); + auto&& read_events_vec = this->read_events(); + auto&& write_events_vec = this->write_events(); + for (std::size_t i = 0; i < write_events_size; ++i) { + read_write_events[i] = read_events_vec[i]; + } + for (std::size_t i = write_events_size, j = 0; i < read_write_size; ++i, ++j) { + read_write_events[i] = write_events_vec[j]; + } + try { + opencl_context.queue().enqueueFillBuffer( + buffer_cl_, static_cast(0), 0, sizeof(T) * this->size(), + &read_write_events, &zero_event); + } catch (const cl::Error& e) { + check_opencl_error("setZero", e); + } + this->add_write_event(zero_event); + } + private: /** * Initializes the OpenCL buffer of this matrix by copying the data from given From 34ad243d2f1113c8abb1c7c014750faa12539b25 Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Tue, 28 Nov 2023 11:43:03 -0500 Subject: [PATCH 06/12] [Jenkins] auto-formatting by clang-format version 10.0.0-4ubuntu1 --- stan/math/opencl/matrix_cl.hpp | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/stan/math/opencl/matrix_cl.hpp b/stan/math/opencl/matrix_cl.hpp index bbb963634c6..54d9129e6ee 100644 --- a/stan/math/opencl/matrix_cl.hpp +++ b/stan/math/opencl/matrix_cl.hpp @@ -519,13 +519,14 @@ class matrix_cl : public matrix_cl_base { for (std::size_t i = 0; i < write_events_size; ++i) { read_write_events[i] = read_events_vec[i]; } - for (std::size_t i = write_events_size, j = 0; i < read_write_size; ++i, ++j) { + for (std::size_t i = write_events_size, j = 0; i < read_write_size; + ++i, ++j) { read_write_events[i] = write_events_vec[j]; } try { - opencl_context.queue().enqueueFillBuffer( - buffer_cl_, static_cast(0), 0, sizeof(T) * this->size(), - &read_write_events, &zero_event); + opencl_context.queue().enqueueFillBuffer(buffer_cl_, static_cast(0), 0, + sizeof(T) * this->size(), + &read_write_events, &zero_event); } catch (const cl::Error& e) { check_opencl_error("setZero", e); } From 7a609a664cac93c1a989824a0c5ae6e84377f6d6 Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Wed, 29 Nov 2023 15:46:19 -0500 Subject: [PATCH 07/12] remove extra setZero() --- stan/math/rev/core/var.hpp | 2 -- 1 file changed, 2 deletions(-) diff --git a/stan/math/rev/core/var.hpp b/stan/math/rev/core/var.hpp index 5553296383e..f6a6d7d7b80 100644 --- a/stan/math/rev/core/var.hpp +++ b/stan/math/rev/core/var.hpp @@ -390,8 +390,6 @@ class var_value> { reverse_pass_callback( [this_vi = this->vi_, other_vi = other.vi_]() mutable { other_vi->adj_ += this_vi->adj_; - // Reset the adjoint for `this` to replicate SoA before assignment - this_vi->adj_.setZero(); }); } From a49035c79a4a30938de484870f158a56a87a9f1d Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Wed, 29 Nov 2023 16:06:29 -0500 Subject: [PATCH 08/12] assign tests y and x should be the same for nullptr case --- test/unit/math/rev/core/var_test.cpp | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/test/unit/math/rev/core/var_test.cpp b/test/unit/math/rev/core/var_test.cpp index 9cb60ba5e80..bad0ecbe9de 100644 --- a/test/unit/math/rev/core/var_test.cpp +++ b/test/unit/math/rev/core/var_test.cpp @@ -935,7 +935,7 @@ TEST_F(AgradRev, assign_nan) { EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); } -TEST_F(AgradRev, assign_nullptr_vari) { +TEST_F(AgradRev, assign_nullptr_var) { using stan::math::var_value; using var_vector = var_value>; using stan::math::var; @@ -954,6 +954,5 @@ TEST_F(AgradRev, assign_nullptr_vari) { x_ans_adj(i) = -(i + 0.1); } EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); - Eigen::VectorXd y_ans_adj = Eigen::VectorXd::Zero(10); - EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); + EXPECT_MATRIX_EQ(x_ans_adj, y.adj()); } From 8449c991370e21ddfe1df1e2e3531510ffcc1261 Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Thu, 30 Nov 2023 15:12:27 -0500 Subject: [PATCH 09/12] change hard_copy to deep_copy, fix async logic in setZero() for matrix_cl, and add test for nan var matrix behaviors --- stan/math/opencl/matrix_cl.hpp | 1 - stan/math/rev/core/arena_matrix.hpp | 2 +- stan/math/rev/core/var.hpp | 14 ++++++++------ test/unit/math/rev/core/var_test.cpp | 28 +++++++++++++++++++++++++++- 4 files changed, 36 insertions(+), 9 deletions(-) diff --git a/stan/math/opencl/matrix_cl.hpp b/stan/math/opencl/matrix_cl.hpp index bbb963634c6..bc71ed09535 100644 --- a/stan/math/opencl/matrix_cl.hpp +++ b/stan/math/opencl/matrix_cl.hpp @@ -509,7 +509,6 @@ class matrix_cl : public matrix_cl_base { return; } cl::Event zero_event; - this->wait_for_read_write_events(); const std::size_t write_events_size = this->write_events().size(); const std::size_t read_events_size = this->read_events().size(); const std::size_t read_write_size = write_events_size + read_events_size; diff --git a/stan/math/rev/core/arena_matrix.hpp b/stan/math/rev/core/arena_matrix.hpp index d6a83f2b6f6..84bb00bca7c 100644 --- a/stan/math/rev/core/arena_matrix.hpp +++ b/stan/math/rev/core/arena_matrix.hpp @@ -134,7 +134,7 @@ class arena_matrix : public Eigen::Map { * @param x the values to write to `this` */ template - void hard_copy(const T& x) { + void deep_copy(const T& x) { Base::operator=(x); } }; diff --git a/stan/math/rev/core/var.hpp b/stan/math/rev/core/var.hpp index 5553296383e..bbac76077f5 100644 --- a/stan/math/rev/core/var.hpp +++ b/stan/math/rev/core/var.hpp @@ -1052,20 +1052,20 @@ class var_value> { return *this; } arena_t> prev_val(vi_->val_.rows(), vi_->val_.cols()); - prev_val.hard_copy(vi_->val_); - vi_->val_.hard_copy(other.val()); + prev_val.deep_copy(vi_->val_); + vi_->val_.deep_copy(other.val()); // no need to change any adjoints - these are just zeros before the reverse // pass reverse_pass_callback( [this_vi = this->vi_, other_vi = other.vi_, prev_val]() mutable { - this_vi->val_.hard_copy(prev_val); + this_vi->val_.deep_copy(prev_val); // we have no way of detecting aliasing between this->vi_->adj_ and // other.vi_->adj_, so we must copy adjoint before reseting to zero // we can reuse prev_val instead of allocating a new matrix - prev_val.hard_copy(this_vi->adj_); + prev_val.deep_copy(this_vi->adj_); this_vi->adj_.setZero(); other_vi->adj_ += prev_val; }); @@ -1074,13 +1074,15 @@ class var_value> { /** * Assignment of another var value, when either both `this` or other does not * contain a plain type. + * @note Here we do not need to use `deep_copy` as the `var_value`'s + * inner `vari_type` holds a view which will call the assignment operator + * that does not perform a placement new. * @tparam S type of the value in the `var_value` to assign * @param other the value to assign * @return this */ template * = nullptr, - require_any_not_plain_type_t* = nullptr, require_not_plain_type_t* = nullptr> inline var_value& operator=(const var_value& other) { // If vi_ is nullptr then the var needs initialized via copy constructor @@ -1105,7 +1107,7 @@ class var_value> { // other.vi_->adj_, so we must copy adjoint before reseting to zero // we can reuse prev_val instead of allocating a new matrix - prev_val.hard_copy(this_vi->adj_); + prev_val = this_vi->adj_; this_vi->adj_.setZero(); other_vi->adj_ += prev_val; }); diff --git a/test/unit/math/rev/core/var_test.cpp b/test/unit/math/rev/core/var_test.cpp index 9cb60ba5e80..3eb8ce77116 100644 --- a/test/unit/math/rev/core/var_test.cpp +++ b/test/unit/math/rev/core/var_test.cpp @@ -911,7 +911,7 @@ TEST_F(AgradRev, matrix_compile_time_conversions) { EXPECT_MATRIX_FLOAT_EQ(x11.val(), rowvec.val()); } -TEST_F(AgradRev, assign_nan) { +TEST_F(AgradRev, assign_nan_varmat) { using stan::math::var_value; using var_vector = var_value>; using stan::math::var; @@ -935,6 +935,32 @@ TEST_F(AgradRev, assign_nan) { EXPECT_MATRIX_EQ(y_ans_adj, y.adj()); } +TEST_F(AgradRev, assign_nan_matvar) { + using stan::math::var; + using var_vector = Eigen::Matrix; + Eigen::VectorXd x_val(10); + for (int i = 0; i < 10; ++i) { + x_val(i) = i + 0.1; + } + var_vector x(x_val); + var_vector y = var_vector(Eigen::Matrix::Constant( + 10, std::numeric_limits::quiet_NaN())); + // need to store y's previous vari pointers + var_vector z = y; + y = stan::math::head(x, 10); + var sigma = 1.0; + var lp = stan::math::normal_lpdf(y, 0, sigma); + lp.grad(); + Eigen::VectorXd x_ans_adj(10); + for (int i = 0; i < 10; ++i) { + x_ans_adj(i) = -(i + 0.1); + } + EXPECT_MATRIX_EQ(x.adj(), x_ans_adj); + Eigen::VectorXd z_ans_adj = Eigen::VectorXd::Zero(10); + EXPECT_MATRIX_EQ(z_ans_adj, z.adj()); +} + + TEST_F(AgradRev, assign_nullptr_vari) { using stan::math::var_value; using var_vector = var_value>; From 42291fb413135948f0338a2bc6fe00823097f2af Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Thu, 30 Nov 2023 16:54:14 -0500 Subject: [PATCH 10/12] add test for matrix_cl.setZero() --- stan/math/opencl/matrix_cl.hpp | 5 +++-- test/unit/math/opencl/matrix_cl_test.cpp | 14 ++++++++++++++ test/unit/math/rev/core/var_test.cpp | 13 +++++++++++++ 3 files changed, 30 insertions(+), 2 deletions(-) diff --git a/stan/math/opencl/matrix_cl.hpp b/stan/math/opencl/matrix_cl.hpp index fb96484003d..90c5e545174 100644 --- a/stan/math/opencl/matrix_cl.hpp +++ b/stan/math/opencl/matrix_cl.hpp @@ -515,10 +515,10 @@ class matrix_cl : public matrix_cl_base { std::vector read_write_events(read_write_size, cl::Event{}); auto&& read_events_vec = this->read_events(); auto&& write_events_vec = this->write_events(); - for (std::size_t i = 0; i < write_events_size; ++i) { + for (std::size_t i = 0; i < read_events_size; ++i) { read_write_events[i] = read_events_vec[i]; } - for (std::size_t i = write_events_size, j = 0; i < read_write_size; + for (std::size_t i = read_events_size, j = 0; j < write_events_size; ++i, ++j) { read_write_events[i] = write_events_vec[j]; } @@ -526,6 +526,7 @@ class matrix_cl : public matrix_cl_base { opencl_context.queue().enqueueFillBuffer(buffer_cl_, static_cast(0), 0, sizeof(T) * this->size(), &read_write_events, &zero_event); + zero_event.wait(); } catch (const cl::Error& e) { check_opencl_error("setZero", e); } diff --git a/test/unit/math/opencl/matrix_cl_test.cpp b/test/unit/math/opencl/matrix_cl_test.cpp index 19e05fda311..92a50033378 100644 --- a/test/unit/math/opencl/matrix_cl_test.cpp +++ b/test/unit/math/opencl/matrix_cl_test.cpp @@ -77,4 +77,18 @@ TEST(MathMatrixCL, assignment) { EXPECT_EQ(nullptr, mat1_cl.buffer()()); } +TEST(MathMatrixCL, setZeroFun) { + using stan::math::matrix_cl; + Eigen::Matrix mat_1; + mat_1 << 1, 2, 3, 4; + matrix_cl mat1_cl(mat_1); + mat1_cl.setZero(); + Eigen::Matrix mat_1_fromcl + = stan::math::from_matrix_cl(mat1_cl); + EXPECT_EQ(mat_1_fromcl(0), 0); + EXPECT_EQ(mat_1_fromcl(1), 0); + EXPECT_EQ(mat_1_fromcl(2), 0); + EXPECT_EQ(mat_1_fromcl(3), 0); +} + #endif diff --git a/test/unit/math/rev/core/var_test.cpp b/test/unit/math/rev/core/var_test.cpp index d291d257ccb..11d7cc2cce6 100644 --- a/test/unit/math/rev/core/var_test.cpp +++ b/test/unit/math/rev/core/var_test.cpp @@ -960,6 +960,19 @@ TEST_F(AgradRev, assign_nan_matvar) { EXPECT_MATRIX_EQ(z_ans_adj, z.adj()); } +/** + * For var and Matrix, we need to make sure + * the tape, when going through reverse mode, leads to the same outcomes. + * In the case where we declare a var without initializing it, aka + * `var_value`, we need to think about what the equivalent + * behavior is for `Eigen::Matrix`. + * When default constructing `Eigen::Matrix` we would have an array of `var` types + * with `nullptr` as the vari. The first assignment to that array would then + * just copy the vari pointer from the other array. This is the behavior we + * want to mimic for `var_value`. So in this test show that + * for uninitialized `var_value`, we can assign it and the + * adjoints are the same as x. + */ TEST_F(AgradRev, assign_nullptr_var) { using stan::math::var_value; using var_vector = var_value>; From 7dcf33b6bbf31df922222d394d9e903ceb4464f8 Mon Sep 17 00:00:00 2001 From: Stan Jenkins Date: Thu, 30 Nov 2023 16:55:09 -0500 Subject: [PATCH 11/12] [Jenkins] auto-formatting by clang-format version 10.0.0-4ubuntu1 --- test/unit/math/rev/core/var_test.cpp | 12 ++++++------ 1 file changed, 6 insertions(+), 6 deletions(-) diff --git a/test/unit/math/rev/core/var_test.cpp b/test/unit/math/rev/core/var_test.cpp index 11d7cc2cce6..025817d408f 100644 --- a/test/unit/math/rev/core/var_test.cpp +++ b/test/unit/math/rev/core/var_test.cpp @@ -966,12 +966,12 @@ TEST_F(AgradRev, assign_nan_matvar) { * In the case where we declare a var without initializing it, aka * `var_value`, we need to think about what the equivalent * behavior is for `Eigen::Matrix`. - * When default constructing `Eigen::Matrix` we would have an array of `var` types - * with `nullptr` as the vari. The first assignment to that array would then - * just copy the vari pointer from the other array. This is the behavior we - * want to mimic for `var_value`. So in this test show that - * for uninitialized `var_value`, we can assign it and the - * adjoints are the same as x. + * When default constructing `Eigen::Matrix` we would have an array + * of `var` types with `nullptr` as the vari. The first assignment to that array + * would then just copy the vari pointer from the other array. This is the + * behavior we want to mimic for `var_value`. So in this test + * show that for uninitialized `var_value`, we can assign it + * and the adjoints are the same as x. */ TEST_F(AgradRev, assign_nullptr_var) { using stan::math::var_value; From bbc701a4b1bf634cccb91da4c41b41af2c3f9638 Mon Sep 17 00:00:00 2001 From: Steve Bronder Date: Thu, 30 Nov 2023 16:55:42 -0500 Subject: [PATCH 12/12] add test for matrix_cl.setZero() --- stan/math/opencl/matrix_cl.hpp | 1 - 1 file changed, 1 deletion(-) diff --git a/stan/math/opencl/matrix_cl.hpp b/stan/math/opencl/matrix_cl.hpp index 90c5e545174..4114c1199c8 100644 --- a/stan/math/opencl/matrix_cl.hpp +++ b/stan/math/opencl/matrix_cl.hpp @@ -526,7 +526,6 @@ class matrix_cl : public matrix_cl_base { opencl_context.queue().enqueueFillBuffer(buffer_cl_, static_cast(0), 0, sizeof(T) * this->size(), &read_write_events, &zero_event); - zero_event.wait(); } catch (const cl::Error& e) { check_opencl_error("setZero", e); }