diff --git a/src/frontend/Typechecker.ml b/src/frontend/Typechecker.ml index 23a05639d..f83371141 100644 --- a/src/frontend/Typechecker.ml +++ b/src/frontend/Typechecker.ml @@ -919,6 +919,8 @@ and check_laplace_fn ~is_cond_dist loc cf tenv id tes = lik_args @ (hbs_arg :: cov_fun_type :: cov_tupl :: control_args) in let return_type = if String.is_suffix id.name ~suffix:"_rng" then UnsizedType.UVector + else if String.is_substring id.name ~substring:"_solve" then + UnsizedType.UTuple [UVector; UMatrix] else UnsizedType.UReal in mk_fun_app ~is_cond_dist ~loc (StanLib (Fun_kind.suffix_from_name id.name)) diff --git a/src/stan_math_signatures/Stan_math_signatures.ml b/src/stan_math_signatures/Stan_math_signatures.ml index 49dedc2d9..c21409acc 100644 --- a/src/stan_math_signatures/Stan_math_signatures.ml +++ b/src/stan_math_signatures/Stan_math_signatures.ml @@ -153,7 +153,8 @@ let is_reduce_sum_fn f = let embedded_laplace_functions = [ (* general fns *) "laplace_marginal"; "laplace_marginal_tol" - ; "laplace_latent_rng"; "laplace_latent_tol_rng"; (* "helpers" *) + ; "laplace_latent_rng"; "laplace_latent_tol_rng"; "laplace_latent_solve" + ; "laplace_latent_solve_tol"; (* "helpers" *) "laplace_marginal_bernoulli_logit_lpmf" ; "laplace_marginal_tol_bernoulli_logit_lpmf" ; "laplace_marginal_neg_binomial_2_log_lpmf" diff --git a/test/integration/bad/embedded_laplace/bad_solve_return.stan b/test/integration/bad/embedded_laplace/bad_solve_return.stan new file mode 100644 index 000000000..2843f1bce --- /dev/null +++ b/test/integration/bad/embedded_laplace/bad_solve_return.stan @@ -0,0 +1,39 @@ +functions { + // specify negative binomial likelihood with mean offset + real ll_function(vector theta, // latent Gaussian + real eta, vector log_ye, // mean offset + array[] int y) { + // observed count + return neg_binomial_2_lpmf(y | exp(log_ye + theta), eta); + } + + // specify covariance function + matrix K_function(array[] vector x, int n_obs, real alpha, real rho) { + matrix[n_obs, n_obs] K = gp_exp_quad_cov(x, alpha, rho); + for (i in 1 : n_obs) + K[i, i] += 1e-8; + return K; + } +} +data { + int n_obs; + int n_coordinates; + array[n_obs] int y; + vector[n_obs] ye; + array[n_obs] vector[n_coordinates] x; +} + +transformed data { + vector[n_obs] log_ye = log(ye); +} +parameters { + real alpha; + real rho; + real eta; +} + +generated quantities { + // laplace_latent_solve returns tuple(vector, matrix), not vector + vector[n_obs] theta = laplace_latent_solve(ll_function, (eta, log_ye, y), 1, + K_function, (x, n_obs, alpha, rho)); +} diff --git a/test/integration/bad/embedded_laplace/bad_solve_tol.stan b/test/integration/bad/embedded_laplace/bad_solve_tol.stan new file mode 100644 index 000000000..fbec3ff86 --- /dev/null +++ b/test/integration/bad/embedded_laplace/bad_solve_tol.stan @@ -0,0 +1,40 @@ +functions { + // specify negative binomial likelihood with mean offset + real ll_function(vector theta, // latent Gaussian + real eta, vector log_ye, // mean offset + array[] int y) { + // observed count + return neg_binomial_2_lpmf(y | exp(log_ye + theta), eta); + } + + // specify covariance function + matrix K_function(array[] vector x, int n_obs, real alpha, real rho) { + matrix[n_obs, n_obs] K = gp_exp_quad_cov(x, alpha, rho); + for (i in 1 : n_obs) + K[i, i] += 1e-8; + return K; + } +} +data { + int n_obs; + int n_coordinates; + array[n_obs] int y; + vector[n_obs] ye; + array[n_obs] vector[n_coordinates] x; +} + +transformed data { + vector[n_obs] log_ye = log(ye); +} +parameters { + real alpha; + real rho; + real eta; +} + +generated quantities { + // _tol variant requires a trailing control-parameter tuple + tuple(vector[n_obs], matrix[n_obs, n_obs]) mean_chol + = laplace_latent_solve_tol(ll_function, (eta, log_ye, y), 1, K_function, + (x, n_obs, alpha, rho)); +} diff --git a/test/integration/bad/embedded_laplace/stanc.expected b/test/integration/bad/embedded_laplace/stanc.expected index 99aa7728e..cf9f91cb6 100644 --- a/test/integration/bad/embedded_laplace/stanc.expected +++ b/test/integration/bad/embedded_laplace/stanc.expected @@ -325,6 +325,35 @@ Semantic error in 'bad_overload.stan', line 2, column 7 to column 51: ------------------------------------------------- Identifier "laplace_marginal_tol_neg_binomial_2_log_lpmf" clashes with a non-overloadable Stan Math library function. +[exit 1] + $ stanc bad_solve_return.stan +Semantic error in 'bad_solve_return.stan', line 37, column 2 to line 38, column 81: + ------------------------------------------------- + 35: generated quantities { + 36: // laplace_latent_solve returns tuple(vector, matrix), not vector + 37: vector[n_obs] theta = laplace_latent_solve(ll_function, (eta, log_ye, y), 1, + ^ + 38: K_function, (x, n_obs, alpha, rho)); + 39: } + ------------------------------------------------- + +Ill-typed assignment statement. +Expected the right hand side to have a type matching the destination (vector). +Instead found type tuple(vector, matrix). +[exit 1] + $ stanc bad_solve_tol.stan +Semantic error in 'bad_solve_tol.stan', line 38, column 8 to line 39, column 56: + ------------------------------------------------- + 36: // _tol variant requires a trailing control-parameter tuple + 37: tuple(vector[n_obs], matrix[n_obs, n_obs]) mean_chol + 38: = laplace_latent_solve_tol(ll_function, (eta, log_ye, y), 1, K_function, + ^ + 39: (x, n_obs, alpha, rho)); + 40: } + ------------------------------------------------- + +Missing control parameter tuple at the end of the call to "laplace_latent_solve_tol". +Expected a tuple of 6 arguments for the control parameters. [exit 1] $ stanc bad_theta0.stan Semantic error in 'bad_theta0.stan', line 43, column 9 to column 16: diff --git a/test/integration/good/code-gen/cpp.expected b/test/integration/good/code-gen/cpp.expected index 4d1a9b293..d94bfa4f3 100644 --- a/test/integration/good/code-gen/cpp.expected +++ b/test/integration/good/code-gen/cpp.expected @@ -7956,13 +7956,15 @@ namespace laplace_functionals_model_namespace { using stan::model::model_base_crtp; using namespace stan::math; stan::math::profile_map profiles__; -static constexpr std::array locations_array__ = +static constexpr std::array locations_array__ = {" (found before start of program)", " (in 'laplace_functionals.stan', line 43, column 2 to column 22)", " (in 'laplace_functionals.stan', line 44, column 2 to column 20)", " (in 'laplace_functionals.stan', line 45, column 2 to column 20)", " (in 'laplace_functionals.stan', line 62, column 2 to line 64, column 67)", " (in 'laplace_functionals.stan', line 66, column 2 to line 72, column 99)", + " (in 'laplace_functionals.stan', line 74, column 2 to line 76, column 65)", + " (in 'laplace_functionals.stan', line 78, column 2 to line 83, column 77)", " (in 'laplace_functionals.stan', line 48, column 2 to column 55)", " (in 'laplace_functionals.stan', line 49, column 2 to column 61)", " (in 'laplace_functionals.stan', line 50, column 2 to column 21)", @@ -7993,6 +7995,12 @@ static constexpr std::array locations_array__ = " (in 'laplace_functionals.stan', line 40, column 2 to column 28)", " (in 'laplace_functionals.stan', line 62, column 9 to column 14)", " (in 'laplace_functionals.stan', line 66, column 9 to column 14)", + " (in 'laplace_functionals.stan', line 74, column 15 to column 20)", + " (in 'laplace_functionals.stan', line 74, column 30 to column 35)", + " (in 'laplace_functionals.stan', line 74, column 37 to column 42)", + " (in 'laplace_functionals.stan', line 78, column 15 to column 20)", + " (in 'laplace_functionals.stan', line 78, column 30 to column 35)", + " (in 'laplace_functionals.stan', line 78, column 37 to column 42)", " (in 'laplace_functionals.stan', line 7, column 4 to column 61)", " (in 'laplace_functionals.stan', line 5, column 34 to line 8, column 3)", " (in 'laplace_functionals.stan', line 12, column 11 to column 16)", @@ -8085,7 +8093,7 @@ ll_function(const T0__& theta_arg__, const T1__& eta, const T2__& // suppress unused var warning (void) DUMMY_VAR__; try { - current_statement__ = 36; + current_statement__ = 44; return stan::math::neg_binomial_2_lpmf(y, stan::math::exp(stan::math::add(log_ye, theta)), eta); } catch (const std::exception& e) { @@ -8117,25 +8125,25 @@ K_function(const T0__& x, const T1__& n_obs, const T2__& alpha, const T3__& // suppress unused var warning (void) DUMMY_VAR__; try { - current_statement__ = 38; + current_statement__ = 46; stan::math::validate_non_negative_index("K", "n_obs", n_obs); - current_statement__ = 39; + current_statement__ = 47; stan::math::validate_non_negative_index("K", "n_obs", n_obs); Eigen::Matrix K = Eigen::Matrix::Constant(n_obs, n_obs, DUMMY_VAR__); - current_statement__ = 40; + current_statement__ = 48; stan::model::assign(K, stan::math::gp_exp_quad_cov(x, alpha, rho), "assigning variable K"); - current_statement__ = 42; + current_statement__ = 50; for (int i = 1; i <= n_obs; ++i) { - current_statement__ = 41; + current_statement__ = 49; stan::model::assign(K, (stan::model::rvalue(K, "K", stan::model::index_uni(i), stan::model::index_uni(i)) + 1e-8), "assigning variable K", stan::model::index_uni(i), stan::model::index_uni(i)); } - current_statement__ = 43; + current_statement__ = 51; return K; } catch (const std::exception& e) { stan::lang::rethrow_located(e, locations_array__[current_statement__]); @@ -8185,29 +8193,29 @@ class laplace_functionals_model final : public model_base_crtp::min(); pos__ = 1; - current_statement__ = 11; + current_statement__ = 13; context__.validate_dims("data initialization", "n_obs", "int", std::vector{}); n_obs = std::numeric_limits::min(); - current_statement__ = 11; + current_statement__ = 13; n_obs = context__.vals_i("n_obs")[(1 - 1)]; - current_statement__ = 12; + current_statement__ = 14; context__.validate_dims("data initialization", "n_coordinates", "int", std::vector{}); n_coordinates = std::numeric_limits::min(); - current_statement__ = 12; + current_statement__ = 14; n_coordinates = context__.vals_i("n_coordinates")[(1 - 1)]; - current_statement__ = 13; + current_statement__ = 15; stan::math::validate_non_negative_index("y", "n_obs", n_obs); - current_statement__ = 14; + current_statement__ = 16; context__.validate_dims("data initialization", "y", "int", std::vector{static_cast(n_obs)}); y = std::vector(n_obs, std::numeric_limits::min()); - current_statement__ = 14; + current_statement__ = 16; y = context__.vals_i("y"); - current_statement__ = 15; + current_statement__ = 17; stan::math::validate_non_negative_index("ye", "n_obs", n_obs); - current_statement__ = 16; + current_statement__ = 18; context__.validate_dims("data initialization", "ye", "double", std::vector{static_cast(n_obs)}); ye_data__ = Eigen::Matrix::Constant(n_obs, @@ -8216,7 +8224,7 @@ class laplace_functionals_model final : public model_base_crtp ye_flat__; - current_statement__ = 16; + current_statement__ = 18; ye_flat__ = context__.vals_r("ye"); pos__ = 1; for (int sym1__ = 1; sym1__ <= n_obs; ++sym1__) { @@ -8225,12 +8233,12 @@ class laplace_functionals_model final : public model_base_crtp{static_cast(n_obs), static_cast(n_coordinates)}); @@ -8239,7 +8247,7 @@ class laplace_functionals_model final : public model_base_crtp::quiet_NaN())); { std::vector x_flat__; - current_statement__ = 19; + current_statement__ = 21; x_flat__ = context__.vals_r("x"); pos__ = 1; for (int sym1__ = 1; sym1__ <= n_coordinates; ++sym1__) { @@ -8251,79 +8259,94 @@ class laplace_functionals_model final : public model_base_crtp{}); rho_location_prior = std::numeric_limits::quiet_NaN(); - current_statement__ = 20; + current_statement__ = 22; rho_location_prior = context__.vals_r("rho_location_prior")[(1 - 1)]; - current_statement__ = 21; + current_statement__ = 23; context__.validate_dims("data initialization", "rho_scale_prior", "double", std::vector{}); rho_scale_prior = std::numeric_limits::quiet_NaN(); - current_statement__ = 21; + current_statement__ = 23; rho_scale_prior = context__.vals_r("rho_scale_prior")[(1 - 1)]; - current_statement__ = 22; + current_statement__ = 24; context__.validate_dims("data initialization", "alpha_location_prior", "double", std::vector{}); alpha_location_prior = std::numeric_limits::quiet_NaN(); - current_statement__ = 22; + current_statement__ = 24; alpha_location_prior = context__.vals_r("alpha_location_prior")[(1 - 1)]; - current_statement__ = 23; + current_statement__ = 25; context__.validate_dims("data initialization", "alpha_scale_prior", "double", std::vector{}); alpha_scale_prior = std::numeric_limits::quiet_NaN(); - current_statement__ = 23; + current_statement__ = 25; alpha_scale_prior = context__.vals_r("alpha_scale_prior")[(1 - 1)]; - current_statement__ = 24; + current_statement__ = 26; stan::math::validate_non_negative_index("log_ye", "n_obs", n_obs); - current_statement__ = 25; + current_statement__ = 27; log_ye_data__ = Eigen::Matrix::Constant(n_obs, std::numeric_limits::quiet_NaN()); new (&log_ye) Eigen::Map>(log_ye_data__.data(), n_obs); - current_statement__ = 25; + current_statement__ = 27; stan::model::assign(log_ye, stan::math::log(ye), "assigning variable log_ye"); - current_statement__ = 26; + current_statement__ = 28; stan::math::validate_non_negative_index("theta_0", "n_obs", n_obs); - current_statement__ = 27; + current_statement__ = 29; theta_0_data__ = Eigen::Matrix::Constant(n_obs, std::numeric_limits::quiet_NaN()); new (&theta_0) Eigen::Map>(theta_0_data__.data(), n_obs); - current_statement__ = 27; + current_statement__ = 29; stan::model::assign(theta_0, stan::math::rep_vector(0.0, n_obs), "assigning variable theta_0"); - current_statement__ = 28; + current_statement__ = 30; tolerance = std::numeric_limits::quiet_NaN(); - current_statement__ = 28; + current_statement__ = 30; tolerance = 1e-6; - current_statement__ = 29; + current_statement__ = 31; max_num_steps = std::numeric_limits::min(); - current_statement__ = 29; + current_statement__ = 31; max_num_steps = 100; - current_statement__ = 30; + current_statement__ = 32; hessian_block_size = std::numeric_limits::min(); - current_statement__ = 30; + current_statement__ = 32; hessian_block_size = 1; - current_statement__ = 31; + current_statement__ = 33; solver = std::numeric_limits::min(); - current_statement__ = 31; + current_statement__ = 33; solver = 1; - current_statement__ = 32; + current_statement__ = 34; max_steps_line_search = std::numeric_limits::min(); - current_statement__ = 32; + current_statement__ = 34; max_steps_line_search = 0; - current_statement__ = 33; + current_statement__ = 35; allow_fallthrough = std::numeric_limits::min(); - current_statement__ = 33; + current_statement__ = 35; allow_fallthrough = 1; - current_statement__ = 34; + current_statement__ = 36; stan::math::validate_non_negative_index("theta", "n_obs", n_obs); - current_statement__ = 35; + current_statement__ = 37; stan::math::validate_non_negative_index("theta2", "n_obs", n_obs); + current_statement__ = 38; + stan::math::validate_non_negative_index("mean_chol.1", "n_obs", n_obs); + current_statement__ = 39; + stan::math::validate_non_negative_index("mean_chol.2", "n_obs", n_obs); + current_statement__ = 40; + stan::math::validate_non_negative_index("mean_chol.2", "n_obs", n_obs); + current_statement__ = 41; + stan::math::validate_non_negative_index("mean_chol_tol.1", "n_obs", + n_obs); + current_statement__ = 42; + stan::math::validate_non_negative_index("mean_chol_tol.2", "n_obs", + n_obs); + current_statement__ = 43; + stan::math::validate_non_negative_index("mean_chol_tol.2", "n_obs", + n_obs); } catch (const std::exception& e) { stan::lang::rethrow_located(e, locations_array__[current_statement__]); } @@ -8373,22 +8396,22 @@ class laplace_functionals_model final : public model_base_crtp(0, lp__); { - current_statement__ = 6; + current_statement__ = 8; lp_accum__.add(stan::math::inv_gamma_lpdf(rho, rho_location_prior, rho_scale_prior)); - current_statement__ = 7; + current_statement__ = 9; lp_accum__.add(stan::math::inv_gamma_lpdf(alpha, alpha_location_prior, alpha_scale_prior)); - current_statement__ = 8; + current_statement__ = 10; lp_accum__.add(stan::math::normal_lpdf(eta, static_cast(0), static_cast(1))); - current_statement__ = 9; + current_statement__ = 11; lp_accum__.add(stan::math::laplace_marginal(ll_function_functor__(), std::forward_as_tuple(eta, log_ye, y), hessian_block_size, K_function_functor__(), std::forward_as_tuple(x, n_obs, alpha, rho), pstream__)); - current_statement__ = 10; + current_statement__ = 12; lp_accum__.add(stan::math::laplace_marginal_tol( ll_function_functor__(), std::forward_as_tuple(eta, log_ye, y), @@ -8441,22 +8464,22 @@ class laplace_functionals_model final : public model_base_crtp(0, lp__); { - current_statement__ = 6; + current_statement__ = 8; lp_accum__.add(stan::math::inv_gamma_lpdf(rho, rho_location_prior, rho_scale_prior)); - current_statement__ = 7; + current_statement__ = 9; lp_accum__.add(stan::math::inv_gamma_lpdf(alpha, alpha_location_prior, alpha_scale_prior)); - current_statement__ = 8; + current_statement__ = 10; lp_accum__.add(stan::math::normal_lpdf(eta, static_cast(0), static_cast(1))); - current_statement__ = 9; + current_statement__ = 11; lp_accum__.add(stan::math::laplace_marginal(ll_function_functor__(), std::forward_as_tuple(eta, log_ye, y), hessian_block_size, K_function_functor__(), std::forward_as_tuple(x, n_obs, alpha, rho), pstream__)); - current_statement__ = 10; + current_statement__ = 12; lp_accum__.add(stan::math::laplace_marginal_tol( ll_function_functor__(), std::forward_as_tuple(eta, log_ye, y), @@ -8552,8 +8575,42 @@ class laplace_functionals_model final : public model_base_crtp, Eigen::Matrix> + mean_chol = + std::tuple, Eigen::Matrix>{ + Eigen::Matrix::Constant(n_obs, + std::numeric_limits::quiet_NaN()), + Eigen::Matrix::Constant(n_obs, n_obs, + std::numeric_limits::quiet_NaN())}; + current_statement__ = 6; + stan::model::assign(mean_chol, + stan::math::laplace_latent_solve(ll_function_functor__(), + std::forward_as_tuple(eta, log_ye, y), hessian_block_size, + K_function_functor__(), + std::forward_as_tuple(x, n_obs, alpha, rho), pstream__), + "assigning variable mean_chol"); + std::tuple, Eigen::Matrix> + mean_chol_tol = + std::tuple, Eigen::Matrix>{ + Eigen::Matrix::Constant(n_obs, + std::numeric_limits::quiet_NaN()), + Eigen::Matrix::Constant(n_obs, n_obs, + std::numeric_limits::quiet_NaN())}; + current_statement__ = 7; + stan::model::assign(mean_chol_tol, + stan::math::laplace_latent_solve_tol(ll_function_functor__(), + std::forward_as_tuple(eta, log_ye, y), hessian_block_size, + K_function_functor__(), + std::forward_as_tuple(x, n_obs, alpha, rho), + std::forward_as_tuple(theta_0, tolerance, max_num_steps, solver, + max_steps_line_search, allow_fallthrough), pstream__), + "assigning variable mean_chol_tol"); out__.write(theta); out__.write(theta2); + out__.write(std::get<0>(mean_chol)); + out__.write(std::get<1>(mean_chol)); + out__.write(std::get<0>(mean_chol_tol)); + out__.write(std::get<1>(mean_chol_tol)); } catch (const std::exception& e) { stan::lang::rethrow_located(e, locations_array__[current_statement__]); } @@ -8635,7 +8692,9 @@ class laplace_functionals_model final : public model_base_crtp{"alpha", "rho", "eta"}; if (emit_transformed_parameters__) {} if (emit_generated_quantities__) { - std::vector temp{"theta", "theta2"}; + std::vector + temp{"theta", "theta2", "mean_chol.1", "mean_chol.2", + "mean_chol_tol.1", "mean_chol_tol.2"}; names__.reserve(names__.size() + temp.size()); names__.insert(names__.end(), temp.begin(), temp.end()); } @@ -8650,7 +8709,13 @@ class laplace_functionals_model final : public model_base_crtp> temp{std::vector{static_cast(n_obs)}, - std::vector{static_cast(n_obs)}}; + std::vector{static_cast(n_obs)}, + std::vector{static_cast(n_obs)}, + std::vector{static_cast(n_obs), + static_cast(n_obs)}, + std::vector{static_cast(n_obs)}, + std::vector{static_cast(n_obs), + static_cast(n_obs)}}; dimss__.reserve(dimss__.size() + temp.size()); dimss__.insert(dimss__.end(), temp.begin(), temp.end()); } @@ -8672,6 +8737,28 @@ class laplace_functionals_model final : public model_base_crtp inline void @@ -8709,7 +8818,7 @@ class laplace_functionals_model final : public model_base_crtp params_i; @@ -8727,7 +8836,7 @@ class laplace_functionals_model final : public model_base_crtp(num_to_write, diff --git a/test/integration/good/code-gen/laplace_functionals.stan b/test/integration/good/code-gen/laplace_functionals.stan index 8948773b2..ba0fa4026 100644 --- a/test/integration/good/code-gen/laplace_functionals.stan +++ b/test/integration/good/code-gen/laplace_functionals.stan @@ -71,4 +71,15 @@ generated quantities { (theta_0, tolerance, max_num_steps, solver, max_steps_line_search, allow_fallthrough)); + tuple(vector[n_obs], matrix[n_obs, n_obs]) mean_chol + = laplace_latent_solve(ll_function, (eta, log_ye, y), hessian_block_size, + K_function, (x, n_obs, alpha, rho)); + + tuple(vector[n_obs], matrix[n_obs, n_obs]) mean_chol_tol + = laplace_latent_solve_tol(ll_function, (eta, log_ye, y), + hessian_block_size, K_function, + (x, n_obs, alpha, rho), + (theta_0, tolerance, max_num_steps, solver, + max_steps_line_search, allow_fallthrough)); + }