add cut_param in rng functions

Author: ygeunkim

inst/include/bvhar/src/bayes/misc/draw.h

@@ -96,12 +96,15 @@ inline void set_grp_id(std::set<int>& own_id, std::set<int> cross_id, const Eige
 // 	}
 // }
 
-inline void cut_param(double& param) {
-	if (param < std::numeric_limits<double>::min() || std::isnan(param)) {
-		param = std::numeric_limits<double>::min();
-	} else if (param > std::numeric_limits<double>::max() || std::isinf(param)) {
-		param = std::numeric_limits<double>::max();
-	}
+inline void cut_param(Eigen::Ref<Eigen::VectorXd> param) {
+	param = param.array().isNaN().select(
+		std::numeric_limits<double>::min(),
+		param.array().isInf().select(
+			std::numeric_limits<double>::max(),
+			param
+		)
+	);
+	param = param.cwiseMax(std::numeric_limits<double>::min()).cwiseMin(std::numeric_limits<double>::max());
 }
 
 // Building Spike-and-slab SD Diagonal Matrix
@@ -784,6 +787,7 @@ inline void dl_local_sparsity(Eigen::VectorXd& local_param, double& dir_concen,
 		cut_param(local_param[i]);
 	}
 	local_param /= local_param.sum();
+	cut_param(local_param);
 }
 
 // Generating Global Parameter of Dirichlet-Laplace Prior

inst/include/bvhar/src/core/common.h

@@ -524,6 +524,14 @@ inline double invgamma_dens(double x, double shp, double scl, bool lg) {
 }
 
 // RNG----------------------------------------
+inline void cut_param(double& param) {
+	if (param < std::numeric_limits<double>::min() || std::isnan(param)) {
+		param = std::numeric_limits<double>::min();
+	} else if (param > std::numeric_limits<double>::max() || std::isinf(param)) {
+		param = std::numeric_limits<double>::max();
+	}
+}
+
 #ifdef USE_RCPP
 inline double bindom_rand(int n, double prob) {
 	return Rf_rbinom(n, prob);
@@ -557,6 +565,7 @@ inline double chisq_rand(double df, BHRNG& rng) {
 }
 
 inline double gamma_rand(double shp, double scl, BHRNG& rng) {
+	cut_param(scl);
 	boost::random::gamma_distribution<> rdist(shp, scl); // 2nd: scale
 	return rdist(rng);
 }

inst/include/bvhar/src/math/random.h

@@ -218,13 +218,17 @@ inline std::vector<Eigen::MatrixXd> sim_mn_iw(const Eigen::MatrixXd& mat_mean, c
 // @param psi Second parameter of GIG
 // @param chi Third parameter of GIG
 inline double sim_gig(double lambda, double psi, double chi, BHRNG& rng) {
+	cut_param(psi);
+	cut_param(chi);
 	boost::random::generalized_inverse_gaussian_distribution<> rdist(lambda, psi, chi);
 	return rdist(rng);
 }
 
 // Generate Inverse Gaussian Distribution
 // This function generates one Inverse Gaussian random number with mu (mean) and lambda (shape).
 inline double sim_invgauss(double mean, double shape, BHRNG& rng) {
+	cut_param(mean);
+	cut_param(shape);
 	boost::random::inverse_gaussian_distribution<> rdist(mean, shape);
 	return rdist(rng);
 }