diff --git a/R/gjamHfunctions.R b/R/gjamHfunctions.R
index 508d24b..8ed14d5 100755
--- a/R/gjamHfunctions.R
+++ b/R/gjamHfunctions.R
@@ -420,7 +420,7 @@ gjamConditionalParameters <- function( output, conditionOn, nsim = 2000 ){
   gindex <- gindex[lkeep,]
   
   wg     <- which(gindex == 1,arr.ind=T)
-  wc     <- matrix(rep(1:M,S*M),S*M,S)[lkeep,]
+  wc     <- matrix(rep(1:M,S*S),S*M,S)[lkeep,]
   rowG   <- wc[wg]
   gindex <- cbind(rowG,wg)
   tmp    <- as.vector( t(outer(colnames(rho)[notOther],
diff --git a/src/cppFns.cpp b/src/cppFns.cpp
index 059fea1..661cdc5 100755
--- a/src/cppFns.cpp
+++ b/src/cppFns.cpp
@@ -115,7 +115,9 @@ arma::mat trMVNmatrixRcpp(arma::mat avec, arma::mat muvec,
       
       cid(0) = cindex;
       idx = idxALLm.col(cindex);
-      sin = arma::inv_sympd(smat.submat(idx, idx));
+      arma::mat smatSub = 0.5 * (smat.submat(idx, idx) + trans(smat.submat(idx, idx)));
+      smatSub = 0.5 * (smatSub + trans(smatSub));
+      sin = arma::inv_sympd(smatSub);
       p1 = trans(smat.submat(idx, cid)) * sin;
       
       m1 = mv[cindex] + dot(p1, (av.elem(idx) - mv.elem(idx)));