diff options
| author | Pjotr Prins | 2018-09-06 12:59:12 +0000 |
|---|---|---|
| committer | Pjotr Prins | 2018-09-06 12:59:12 +0000 |
| commit | b40601c6e65b10e070dc3d025ffa50850b7d4ebb (patch) | |
| tree | 9cd3f7138af5a2cb077c03edf66aa4638093b4d1 /src/gemma.cpp | |
| parent | cdf407bd7994dbe41a952bf29cacc1a2ca9c722e (diff) | |
| download | pangemma-b40601c6e65b10e070dc3d025ffa50850b7d4ebb.tar.gz | |
Sometimes a value gets negative zero causing a NaN on the sqrt.
Closes #61
Diffstat (limited to 'src/gemma.cpp')
| -rw-r--r-- | src/gemma.cpp | 12 |
1 files changed, 6 insertions, 6 deletions
diff --git a/src/gemma.cpp b/src/gemma.cpp index cb01ee0..f928977 100644 --- a/src/gemma.cpp +++ b/src/gemma.cpp @@ -2595,7 +2595,7 @@ void GEMMA::BatchRun(PARAM &cPar) { if (wi <= 0 || wj <= 0) { d = 0; } else { - d /= sqrt(wi * wj); + d /= safe_sqrt(wi * wj); } gsl_matrix_set(G, i, j, d); if (j != i) { @@ -2623,7 +2623,7 @@ void GEMMA::BatchRun(PARAM &cPar) { if (wi <= 0) { wi = 0; } else { - wi = sqrt(wi); + wi = safe_sqrt(wi); } gsl_vector_view Urow = gsl_matrix_row(U, i); gsl_vector_scale(&Urow.vector, wi); @@ -3418,7 +3418,7 @@ void GEMMA::WriteLog(int argc, char **argv, PARAM &cPar) { for (size_t j = 0; j <= i; j++) { c = (2 * cPar.n_ph - min(i, j) + 1) * min(i, j) / 2 + max(i, j) - min(i, j); - outfile << tab(j) << sqrt(cPar.VVg_remle_null[c]); + outfile << tab(j) << safe_sqrt(cPar.VVg_remle_null[c]); } outfile << endl; } @@ -3436,7 +3436,7 @@ void GEMMA::WriteLog(int argc, char **argv, PARAM &cPar) { for (size_t j = 0; j <= i; j++) { c = (2 * cPar.n_ph - min(i, j) + 1) * min(i, j) / 2 + max(i, j) - min(i, j); - outfile << tab(j) << sqrt(cPar.VVe_remle_null[c]); + outfile << tab(j) << safe_sqrt(cPar.VVe_remle_null[c]); } outfile << endl; } @@ -3455,7 +3455,7 @@ void GEMMA::WriteLog(int argc, char **argv, PARAM &cPar) { for (size_t j = 0; j <= i; j++) { c = (2 * cPar.n_ph - min(i, j) + 1) * min(i, j) / 2 + max(i, j) - min(i, j); - outfile << tab(j) << sqrt(cPar.VVg_mle_null[c]); + outfile << tab(j) << safe_sqrt(cPar.VVg_mle_null[c]); } outfile << endl; } @@ -3473,7 +3473,7 @@ void GEMMA::WriteLog(int argc, char **argv, PARAM &cPar) { for (size_t j = 0; j <= i; j++) { c = (2 * cPar.n_ph - min(i, j) + 1) * min(i, j) / 2 + max(i, j) - min(i, j); - outfile << tab(j) << sqrt(cPar.VVe_mle_null[c]); + outfile << tab(j) << safe_sqrt(cPar.VVe_mle_null[c]); } outfile << endl; } |