1 #ifdef STAN_OPENCL
2 #include <stan/math/prim.hpp>
3 #include <test/unit/math/opencl/util.hpp>
4 #include <test/unit/util.hpp>
5 #include <gtest/gtest.h>
6 #include <algorithm>
7 
TEST(MathMatrixOpenCLPrim,non_matching_dim_excpetion)8 TEST(MathMatrixOpenCLPrim, non_matching_dim_excpetion) {
9   stan::math::matrix_d m0(5, 3);
10   stan::math::matrix_d m1(2, 6);
11 
12   stan::math::matrix_cl<double> m0_cl(m0);
13   stan::math::matrix_cl<double> m1_cl(m1);
14   EXPECT_THROW(m0_cl * m1_cl, std::invalid_argument);
15   EXPECT_THROW(stan::math::multiply(m0_cl, m1_cl), std::invalid_argument);
16 }
17 
TEST(MathMatrixOpenCLPrim,matrix_vector_tri_small)18 TEST(MathMatrixOpenCLPrim, matrix_vector_tri_small) {
19   auto m = stan::math::matrix_d::Random(4, 5).eval();
20   auto v = stan::math::vector_d::Random(5).eval();
21   stan::math::matrix_d m0(4, 1);
22   stan::math::matrix_d m0_cl_res(4, 1);
23 
24   stan::math::matrix_cl<double> v_cl(v);
25   stan::math::matrix_cl<double> m_cl(m);
26   stan::math::matrix_cl<double> m0_cl(4, 1);
27 
28   m0 = m * v.triangularView<Eigen::Lower>();
29   m_cl.view(stan::math::matrix_cl_view::Entire);
30   v_cl.view(stan::math::matrix_cl_view::Lower);
31   m0_cl = m_cl * v_cl;
32   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
33   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
34 
35   m0 = m.triangularView<Eigen::Lower>() * v;
36   m_cl.view(stan::math::matrix_cl_view::Lower);
37   v_cl.view(stan::math::matrix_cl_view::Entire);
38   m0_cl = m_cl * v_cl;
39   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
40   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
41 
42   m0 = m * v.triangularView<Eigen::Upper>();
43   m_cl.view(stan::math::matrix_cl_view::Entire);
44   v_cl.view(stan::math::matrix_cl_view::Upper);
45   m0_cl = m_cl * v_cl;
46   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
47   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
48 
49   m0 = m.triangularView<Eigen::Upper>() * v;
50   m_cl.view(stan::math::matrix_cl_view::Upper);
51   v_cl.view(stan::math::matrix_cl_view::Entire);
52   m0_cl = m_cl * v_cl;
53   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
54   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
55 
56   // the cast is needed because operator* for for two triangular views does not
57   // exist
58   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Lower>())
59        * v.triangularView<Eigen::Lower>();
60   m_cl.view(stan::math::matrix_cl_view::Lower);
61   v_cl.view(stan::math::matrix_cl_view::Lower);
62   m0_cl = m_cl * v_cl;
63   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
64   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
65 
66   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Upper>())
67        * v.triangularView<Eigen::Upper>();
68   m_cl.view(stan::math::matrix_cl_view::Upper);
69   v_cl.view(stan::math::matrix_cl_view::Upper);
70   m0_cl = m_cl * v_cl;
71   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
72   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
73 
74   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Lower>())
75        * v.triangularView<Eigen::Upper>();
76   m_cl.view(stan::math::matrix_cl_view::Lower);
77   v_cl.view(stan::math::matrix_cl_view::Upper);
78   m0_cl = m_cl * v_cl;
79   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
80   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
81 
82   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Upper>())
83        * v.triangularView<Eigen::Lower>();
84   m_cl.view(stan::math::matrix_cl_view::Upper);
85   v_cl.view(stan::math::matrix_cl_view::Lower);
86   m0_cl = m_cl * v_cl;
87   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
88   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
89 }
90 
TEST(MathMatrixOpenCLPrim,matrix_vector_tri_big)91 TEST(MathMatrixOpenCLPrim, matrix_vector_tri_big) {
92   auto m = stan::math::matrix_d::Random(400, 600).eval();
93   auto v = stan::math::vector_d::Random(600).eval();
94   stan::math::matrix_d m0(400, 1);
95   stan::math::matrix_d m0_cl_res(400, 1);
96 
97   stan::math::matrix_cl<double> v_cl(v);
98   stan::math::matrix_cl<double> m_cl(m);
99   stan::math::matrix_cl<double> m0_cl(400, 1);
100 
101   m0 = m * v.triangularView<Eigen::Lower>();
102   m_cl.view(stan::math::matrix_cl_view::Entire);
103   v_cl.view(stan::math::matrix_cl_view::Lower);
104   m0_cl = m_cl * v_cl;
105   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
106   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
107 
108   m0 = m.triangularView<Eigen::Lower>() * v;
109   m_cl.view(stan::math::matrix_cl_view::Lower);
110   v_cl.view(stan::math::matrix_cl_view::Entire);
111   m0_cl = m_cl * v_cl;
112   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
113   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
114 
115   m0 = m * v.triangularView<Eigen::Upper>();
116   m_cl.view(stan::math::matrix_cl_view::Entire);
117   v_cl.view(stan::math::matrix_cl_view::Upper);
118   m0_cl = m_cl * v_cl;
119   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
120   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
121 
122   m0 = m.triangularView<Eigen::Upper>() * v;
123   m_cl.view(stan::math::matrix_cl_view::Upper);
124   v_cl.view(stan::math::matrix_cl_view::Entire);
125   m0_cl = m_cl * v_cl;
126   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
127   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
128 
129   // the cast is needed because operator* for for two triangular views does not
130   // exist
131   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Lower>())
132        * v.triangularView<Eigen::Lower>();
133   m_cl.view(stan::math::matrix_cl_view::Lower);
134   v_cl.view(stan::math::matrix_cl_view::Lower);
135   m0_cl = m_cl * v_cl;
136   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
137   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
138 
139   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Upper>())
140        * v.triangularView<Eigen::Upper>();
141   m_cl.view(stan::math::matrix_cl_view::Upper);
142   v_cl.view(stan::math::matrix_cl_view::Upper);
143   m0_cl = m_cl * v_cl;
144   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
145   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
146 
147   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Lower>())
148        * v.triangularView<Eigen::Upper>();
149   m_cl.view(stan::math::matrix_cl_view::Lower);
150   v_cl.view(stan::math::matrix_cl_view::Upper);
151   m0_cl = m_cl * v_cl;
152   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
153   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
154 
155   m0 = static_cast<Eigen::MatrixXd>(m.triangularView<Eigen::Upper>())
156        * v.triangularView<Eigen::Lower>();
157   m_cl.view(stan::math::matrix_cl_view::Upper);
158   v_cl.view(stan::math::matrix_cl_view::Lower);
159   m0_cl = m_cl * v_cl;
160   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
161   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
162 }
163 
TEST(MathMatrixOpenCLPrim,row_vector_matrix_tri_small)164 TEST(MathMatrixOpenCLPrim, row_vector_matrix_tri_small) {
165   auto m = stan::math::matrix_d::Random(5, 4).eval();
166   auto rv = stan::math::row_vector_d::Random(5).eval();
167   stan::math::matrix_d m0(1, 4);
168   stan::math::matrix_d m0_cl_res(1, 4);
169 
170   stan::math::matrix_cl<double> m_cl(m);
171   stan::math::matrix_cl<double> rv_cl(rv);
172   stan::math::matrix_cl<double> m0_cl(1, 4);
173 
174   m0 = rv * m;
175   m0_cl = rv_cl * m_cl;
176   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
177   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
178 
179   m0 = rv * m.triangularView<Eigen::Lower>();
180   rv_cl.view(stan::math::matrix_cl_view::Entire);
181   m_cl.view(stan::math::matrix_cl_view::Lower);
182   m0_cl = rv_cl * m_cl;
183   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
184   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
185 
186   m0 = rv.triangularView<Eigen::Lower>() * m;
187   rv_cl.view(stan::math::matrix_cl_view::Lower);
188   m_cl.view(stan::math::matrix_cl_view::Entire);
189   m0_cl = rv_cl * m_cl;
190   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
191   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
192 
193   m0 = rv * m.triangularView<Eigen::Upper>();
194   rv_cl.view(stan::math::matrix_cl_view::Entire);
195   m_cl.view(stan::math::matrix_cl_view::Upper);
196   m0_cl = rv_cl * m_cl;
197   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
198   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
199 
200   m0 = rv.triangularView<Eigen::Upper>() * m;
201   rv_cl.view(stan::math::matrix_cl_view::Upper);
202   m_cl.view(stan::math::matrix_cl_view::Entire);
203   m0_cl = rv_cl * m_cl;
204   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
205   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
206 
207   // the cast is needed because operator* for for two triangular views does not
208   // exist
209   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Lower>())
210        * m.triangularView<Eigen::Lower>();
211   rv_cl.view(stan::math::matrix_cl_view::Lower);
212   m_cl.view(stan::math::matrix_cl_view::Lower);
213   m0_cl = rv_cl * m_cl;
214   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
215   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
216 
217   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Upper>())
218        * m.triangularView<Eigen::Upper>();
219   rv_cl.view(stan::math::matrix_cl_view::Upper);
220   m_cl.view(stan::math::matrix_cl_view::Upper);
221   m0_cl = rv_cl * m_cl;
222   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
223   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
224 
225   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Lower>())
226        * m.triangularView<Eigen::Upper>();
227   rv_cl.view(stan::math::matrix_cl_view::Lower);
228   m_cl.view(stan::math::matrix_cl_view::Upper);
229   m0_cl = rv_cl * m_cl;
230   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
231   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
232 
233   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Upper>())
234        * m.triangularView<Eigen::Lower>();
235   rv_cl.view(stan::math::matrix_cl_view::Upper);
236   m_cl.view(stan::math::matrix_cl_view::Lower);
237   m0_cl = rv_cl * m_cl;
238   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
239   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
240 }
241 
TEST(MathMatrixOpenCLPrim,row_vector_matrix_tri_big)242 TEST(MathMatrixOpenCLPrim, row_vector_matrix_tri_big) {
243   auto m = stan::math::matrix_d::Random(600, 400).eval();
244   auto rv = stan::math::row_vector_d::Random(600).eval();
245   stan::math::matrix_d m0(1, 400);
246   stan::math::matrix_d m0_cl_res(1, 400);
247 
248   stan::math::matrix_cl<double> m_cl(m);
249   stan::math::matrix_cl<double> rv_cl(rv);
250   stan::math::matrix_cl<double> m0_cl(1, 400);
251 
252   m0 = rv * m;
253   m0_cl = rv_cl * m_cl;
254   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
255   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
256 
257   m0 = rv * m.triangularView<Eigen::Lower>();
258   rv_cl.view(stan::math::matrix_cl_view::Entire);
259   m_cl.view(stan::math::matrix_cl_view::Lower);
260   m0_cl = rv_cl * m_cl;
261   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
262   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
263 
264   m0 = rv.triangularView<Eigen::Lower>() * m;
265   rv_cl.view(stan::math::matrix_cl_view::Lower);
266   m_cl.view(stan::math::matrix_cl_view::Entire);
267   m0_cl = rv_cl * m_cl;
268   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
269   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
270 
271   m0 = rv * m.triangularView<Eigen::Upper>();
272   rv_cl.view(stan::math::matrix_cl_view::Entire);
273   m_cl.view(stan::math::matrix_cl_view::Upper);
274   m0_cl = rv_cl * m_cl;
275   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
276   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
277 
278   m0 = rv.triangularView<Eigen::Upper>() * m;
279   rv_cl.view(stan::math::matrix_cl_view::Upper);
280   m_cl.view(stan::math::matrix_cl_view::Entire);
281   m0_cl = rv_cl * m_cl;
282   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
283   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
284 
285   // the cast is needed because operator* for for two triangular views does not
286   // exist
287   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Lower>())
288        * m.triangularView<Eigen::Lower>();
289   rv_cl.view(stan::math::matrix_cl_view::Lower);
290   m_cl.view(stan::math::matrix_cl_view::Lower);
291   m0_cl = rv_cl * m_cl;
292   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
293   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
294 
295   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Upper>())
296        * m.triangularView<Eigen::Upper>();
297   rv_cl.view(stan::math::matrix_cl_view::Upper);
298   m_cl.view(stan::math::matrix_cl_view::Upper);
299   m0_cl = rv_cl * m_cl;
300   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
301   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
302 
303   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Lower>())
304        * m.triangularView<Eigen::Upper>();
305   rv_cl.view(stan::math::matrix_cl_view::Lower);
306   m_cl.view(stan::math::matrix_cl_view::Upper);
307   m0_cl = rv_cl * m_cl;
308   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
309   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
310 
311   m0 = static_cast<Eigen::MatrixXd>(rv.triangularView<Eigen::Upper>())
312        * m.triangularView<Eigen::Lower>();
313   rv_cl.view(stan::math::matrix_cl_view::Upper);
314   m_cl.view(stan::math::matrix_cl_view::Lower);
315   m0_cl = rv_cl * m_cl;
316   m0_cl_res = stan::math::from_matrix_cl(m0_cl);
317   EXPECT_MATRIX_NEAR(m0, m0_cl_res, 1e-10);
318 }
319 
TEST(MathMatrixOpenCLPrim,multiply_small)320 TEST(MathMatrixOpenCLPrim, multiply_small) {
321   using stan::math::multiply;
322   auto m1 = stan::math::matrix_d::Random(3, 3).eval();
323   auto m2 = stan::math::matrix_d::Random(3, 3).eval();
324   stan::math::matrix_d m3_cl_res(3, 3);
325 
326   stan::math::matrix_cl<double> m11(m1);
327   stan::math::matrix_cl<double> m22(m2);
328 
329   auto m3 = (m1 * m2).eval();
330 
331   auto m33 = m11 * m22;
332 
333   m3_cl_res = stan::math::from_matrix_cl(m33);
334 
335   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
336 }
337 
TEST(MathMatrixOpenCLPrim,multiply_big)338 TEST(MathMatrixOpenCLPrim, multiply_big) {
339   using stan::math::multiply;
340   int size = 512;
341   auto m1 = stan::math::matrix_d::Random(size, size).eval();
342   auto m2 = stan::math::matrix_d::Random(size, size).eval();
343   stan::math::matrix_d m3_cl_res(size, size);
344 
345   stan::math::matrix_cl<double> m11(m1);
346   stan::math::matrix_cl<double> m22(m2);
347 
348   auto m3 = (m1 * m2).eval();
349 
350   auto m33 = m11 * m22;
351 
352   m3_cl_res = stan::math::from_matrix_cl(m33);
353 
354   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
355 }
356 
TEST(MathMatrixOpenCLPrim,lower_tri_rect_multiply_small)357 TEST(MathMatrixOpenCLPrim, lower_tri_rect_multiply_small) {
358   auto m1 = stan::math::matrix_d::Random(3, 3).eval();
359   auto m2 = stan::math::matrix_d::Random(3, 3).eval();
360   stan::math::matrix_d m3_cl_res(3, 3);
361 
362   stan::math::matrix_cl<double> m11(m1, stan::math::matrix_cl_view::Lower);
363   stan::math::matrix_cl<double> m22(m2);
364 
365   auto m3 = (m1.triangularView<Eigen::Lower>() * m2).eval();
366 
367   auto m33 = m11 * m22;
368 
369   m3_cl_res = stan::math::from_matrix_cl(m33);
370 
371   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
372 }
373 
TEST(MathMatrixOpenCLPrim,lower_tri_rect_multiply_big)374 TEST(MathMatrixOpenCLPrim, lower_tri_rect_multiply_big) {
375   int size = 512;
376   auto m1 = stan::math::matrix_d::Random(size, size).eval();
377   auto m2 = stan::math::matrix_d::Random(size, size).eval();
378   stan::math::matrix_d m3_cl_res(size, size);
379 
380   stan::math::matrix_cl<double> m11(m1, stan::math::matrix_cl_view::Lower);
381   stan::math::matrix_cl<double> m22(m2);
382 
383   auto m3 = (m1.triangularView<Eigen::Lower>() * m2).eval();
384 
385   auto m33 = m11 * m22;
386 
387   m3_cl_res = stan::math::from_matrix_cl(m33);
388 
389   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
390 }
391 
TEST(MathMatrixOpenCLPrim,lower_tri_rect_multiply_big_rect)392 TEST(MathMatrixOpenCLPrim, lower_tri_rect_multiply_big_rect) {
393   int size = 321;
394   auto m1 = stan::math::matrix_d::Random(size, size).eval();
395   auto m2 = stan::math::matrix_d::Random(size, size * 3).eval();
396   stan::math::matrix_d m3_cl_res(size, size * 3);
397 
398   stan::math::matrix_cl<double> m11(m1, stan::math::matrix_cl_view::Lower);
399   stan::math::matrix_cl<double> m22(m2);
400 
401   auto m3 = (m1.triangularView<Eigen::Lower>() * m2).eval();
402 
403   auto m33 = m11 * m22;
404 
405   m3_cl_res = stan::math::from_matrix_cl(m33);
406 
407   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
408 }
409 
TEST(MathMatrixOpenCLPrim,upper_tri_rect_multiply_small)410 TEST(MathMatrixOpenCLPrim, upper_tri_rect_multiply_small) {
411   auto m1 = stan::math::matrix_d::Random(3, 3).eval();
412   auto m2 = stan::math::matrix_d::Random(3, 3).eval();
413   stan::math::matrix_d m3_cl_res(3, 3);
414 
415   stan::math::matrix_cl<double> m11(m1, stan::math::matrix_cl_view::Upper);
416   stan::math::matrix_cl<double> m22(m2);
417 
418   auto m3 = (m1.triangularView<Eigen::Upper>() * m2).eval();
419 
420   auto m33 = m11 * m22;
421 
422   m3_cl_res = stan::math::from_matrix_cl(m33);
423 
424   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
425 }
426 
TEST(MathMatrixOpenCLPrim,upper_tri_rect_multiply_big)427 TEST(MathMatrixOpenCLPrim, upper_tri_rect_multiply_big) {
428   int size = 472;
429   auto m1 = stan::math::matrix_d::Random(size, size).eval();
430   auto m2 = stan::math::matrix_d::Random(size, size).eval();
431   stan::math::matrix_d m3_cl_res(size, size);
432 
433   stan::math::matrix_cl<double> m11(m1, stan::math::matrix_cl_view::Upper);
434   stan::math::matrix_cl<double> m22(m2);
435 
436   auto m3 = (m1.triangularView<Eigen::Upper>() * m2).eval();
437 
438   auto m33 = m11 * m22;
439 
440   m3_cl_res = stan::math::from_matrix_cl(m33);
441 
442   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
443 }
444 
TEST(MathMatrixOpenCLPrim,upper_tri_rect_multiply_big_rect)445 TEST(MathMatrixOpenCLPrim, upper_tri_rect_multiply_big_rect) {
446   int size = 463;
447   auto m1 = stan::math::matrix_d::Random(size, size).eval();
448   auto m2 = stan::math::matrix_d::Random(size, size * 3).eval();
449   stan::math::matrix_d m3_cl_res(size, size * 3);
450 
451   stan::math::matrix_cl<double> m11(m1, stan::math::matrix_cl_view::Upper);
452   stan::math::matrix_cl<double> m22(m2);
453 
454   auto m3 = (m1.triangularView<Eigen::Upper>() * m2).eval();
455 
456   auto m33 = m11 * m22;
457 
458   m3_cl_res = stan::math::from_matrix_cl(m33);
459 
460   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
461 }
462 
TEST(MathMatrixOpenCLPrim,rect_lower_tri_multiply_small)463 TEST(MathMatrixOpenCLPrim, rect_lower_tri_multiply_small) {
464   auto m1 = stan::math::matrix_d::Random(3, 3).eval();
465   auto m2 = stan::math::matrix_d::Random(3, 3).eval();
466   stan::math::matrix_d m3_cl_res(3, 3);
467 
468   stan::math::matrix_cl<double> m11(m1);
469   stan::math::matrix_cl<double> m22(m2, stan::math::matrix_cl_view::Lower);
470 
471   auto m3 = (m1 * m2.triangularView<Eigen::Lower>()).eval();
472 
473   auto m33 = m11 * m22;
474 
475   m3_cl_res = stan::math::from_matrix_cl(m33);
476 
477   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
478 }
479 
TEST(MathMatrixOpenCLPrim,rect_lower_tri_multiply_big)480 TEST(MathMatrixOpenCLPrim, rect_lower_tri_multiply_big) {
481   int size = 451;
482   auto m1 = stan::math::matrix_d::Random(size, size).eval();
483   auto m2 = stan::math::matrix_d::Random(size, size).eval();
484   stan::math::matrix_d m3_cl_res(size, size);
485 
486   stan::math::matrix_cl<double> m11(m1);
487   stan::math::matrix_cl<double> m22(m2, stan::math::matrix_cl_view::Lower);
488 
489   auto m3 = (m1 * m2.triangularView<Eigen::Lower>()).eval();
490 
491   auto m33 = m11 * m22;
492 
493   m3_cl_res = stan::math::from_matrix_cl(m33);
494 
495   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
496 }
497 
TEST(MathMatrixOpenCLPrim,rect_lower_tri_multiply_big_rect)498 TEST(MathMatrixOpenCLPrim, rect_lower_tri_multiply_big_rect) {
499   int size = 444;
500   auto m1 = stan::math::matrix_d::Random(size, size).eval();
501   auto m2 = stan::math::matrix_d::Random(size, size * 3).eval();
502   stan::math::matrix_d m3_cl_res(size, size * 3);
503 
504   stan::math::matrix_cl<double> m11(m1);
505   stan::math::matrix_cl<double> m22(m2, stan::math::matrix_cl_view::Lower);
506 
507   auto m3 = (m1 * m2.triangularView<Eigen::Lower>()).eval();
508 
509   auto m33 = m11 * m22;
510 
511   m3_cl_res = stan::math::from_matrix_cl(m33);
512 
513   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
514 }
515 
TEST(MathMatrixOpenCLPrim,rect_upper_tri_multiply_small)516 TEST(MathMatrixOpenCLPrim, rect_upper_tri_multiply_small) {
517   auto m1 = stan::math::matrix_d::Random(3, 3).eval();
518   auto m2 = stan::math::matrix_d::Random(3, 3).eval();
519   stan::math::matrix_d m3_cl_res(3, 3);
520 
521   stan::math::matrix_cl<double> m11(m1);
522   stan::math::matrix_cl<double> m22(m2, stan::math::matrix_cl_view::Upper);
523 
524   auto m3 = (m1 * m2.triangularView<Eigen::Upper>()).eval();
525 
526   auto m33 = m11 * m22;
527 
528   m3_cl_res = stan::math::from_matrix_cl(m33);
529 
530   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
531 }
532 
TEST(MathMatrixOpenCLPrim,rect_upper_tri_multiply_big)533 TEST(MathMatrixOpenCLPrim, rect_upper_tri_multiply_big) {
534   int size = 468;
535   auto m1 = stan::math::matrix_d::Random(size, size).eval();
536   auto m2 = stan::math::matrix_d::Random(size, size).eval();
537   stan::math::matrix_d m3_cl_res(size, size);
538 
539   stan::math::matrix_cl<double> m11(m1);
540   stan::math::matrix_cl<double> m22(m2, stan::math::matrix_cl_view::Upper);
541 
542   auto m3 = (m1 * m2.triangularView<Eigen::Upper>()).eval();
543 
544   auto m33 = m11 * m22;
545 
546   m3_cl_res = stan::math::from_matrix_cl(m33);
547 
548   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
549 }
550 
TEST(MathMatrixOpenCLPrim,rect_upper_tri_multiply_big_rect)551 TEST(MathMatrixOpenCLPrim, rect_upper_tri_multiply_big_rect) {
552   int size = 345;
553   auto m1 = stan::math::matrix_d::Random(size * 3, size).eval();
554   auto m2 = stan::math::matrix_d::Random(size, size).eval();
555   stan::math::matrix_d m3_cl_res(size * 3, size);
556 
557   stan::math::matrix_cl<double> m11(m1);
558   stan::math::matrix_cl<double> m22(m2, stan::math::matrix_cl_view::Upper);
559 
560   auto m3 = (m1 * m2.triangularView<Eigen::Upper>()).eval();
561 
562   auto m33 = m11 * m22;
563 
564   m3_cl_res = stan::math::from_matrix_cl(m33);
565 
566   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
567 }
568 
TEST(MathMatrixOpenCLPrim,multiply_big_split_4)569 TEST(MathMatrixOpenCLPrim, multiply_big_split_4) {
570   using stan::math::multiply;
571   int size = 512;
572   auto m1 = stan::math::matrix_d::Random(size, size * 2).eval();
573   auto m2 = stan::math::matrix_d::Random(size * 2, size).eval();
574   stan::math::matrix_d m3_cl_res(size, size);
575 
576   stan::math::matrix_cl<double> m11(m1);
577   stan::math::matrix_cl<double> m22(m2);
578 
579   auto m3 = (m1 * m2).eval();
580 
581   auto m33 = stan::math::multiply(m11, m22);
582 
583   m3_cl_res = stan::math::from_matrix_cl(m33);
584 
585   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
586 }
587 
TEST(MathMatrixOpenCLPrim,multiply_big_split_11)588 TEST(MathMatrixOpenCLPrim, multiply_big_split_11) {
589   using stan::math::multiply;
590   int size = 433;
591   auto m1 = stan::math::matrix_d::Random(size, size * 11).eval();
592   auto m2 = stan::math::matrix_d::Random(size * 11, size).eval();
593   stan::math::matrix_d m3_cl_res(size, size);
594 
595   stan::math::matrix_cl<double> m11(m1);
596   stan::math::matrix_cl<double> m22(m2);
597 
598   auto m3 = (m1 * m2).eval();
599 
600   auto m33 = stan::math::multiply(m11, m22);
601 
602   m3_cl_res = stan::math::from_matrix_cl(m33);
603 
604   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
605 }
606 
TEST(MathMatrixOpenCLPrim,multiply_small_split_big)607 TEST(MathMatrixOpenCLPrim, multiply_small_split_big) {
608   using stan::math::multiply;
609   int size = 32;
610   auto m1 = stan::math::matrix_d::Random(size, size * 10).eval();
611   auto m2 = stan::math::matrix_d::Random(size * 10, size).eval();
612   stan::math::matrix_d m3_cl_res(size, size);
613 
614   stan::math::matrix_cl<double> m11(m1);
615   stan::math::matrix_cl<double> m22(m2);
616 
617   auto m3 = (m1 * m2).eval();
618 
619   auto m33 = stan::math::multiply(m11, m22);
620 
621   m3_cl_res = stan::math::from_matrix_cl(m33);
622 
623   EXPECT_MATRIX_NEAR(m3, m3_cl_res, 1e-10);
624 }
625 #endif
626