diff --git a/lib/linalg/README b/lib/linalg/README
index 6c04225d17..61776e6653 100644
--- a/lib/linalg/README
+++ b/lib/linalg/README
@@ -11,7 +11,7 @@ resulting library will follow the Fortran binary conventions.
 
 Note that this is an *incomplete* subset of full BLAS/LAPACK.
 
-The files correspond to LAPACK version 3.11.0.
+The files correspond to LAPACK version 3.12.0.
 
 You should only need to build and use the library in this directory if
 you want to build LAMMPS with one of the listed packages AND you do not
diff --git a/lib/linalg/dbdsqr.cpp b/lib/linalg/dbdsqr.cpp
index 59498b4ae6..988223d50e 100644
--- a/lib/linalg/dbdsqr.cpp
+++ b/lib/linalg/dbdsqr.cpp
@@ -40,7 +40,7 @@ int dbdsqr_(char *uplo, integer *n, integer *ncvt, integer *nru, integer *ncc, d
     integer oldll;
     doublereal shift, sigmn, oldsn;
     extern int dswap_(integer *, doublereal *, integer *, doublereal *, integer *);
-    doublereal sminl, sigmx;
+    doublereal sigmx;
     logical lower;
     extern int dlasq1_(integer *, doublereal *, doublereal *, doublereal *, integer *),
         dlasv2_(doublereal *, doublereal *, doublereal *, doublereal *, doublereal *, doublereal *,
@@ -141,7 +141,7 @@ int dbdsqr_(char *uplo, integer *n, integer *ncvt, integer *nru, integer *ncc, d
         d__2 = smax, d__3 = (d__1 = e[i__], abs(d__1));
         smax = max(d__2, d__3);
     }
-    sminl = 0.;
+    smin = 0.;
     if (tol >= 0.) {
         sminoa = abs(d__[1]);
         if (sminoa == 0.) {
@@ -185,7 +185,6 @@ L60:
         d__[m] = 0.;
     }
     smax = (d__1 = d__[m], abs(d__1));
-    smin = smax;
     i__1 = m - 1;
     for (lll = 1; lll <= i__1; ++lll) {
         ll = m - lll;
@@ -197,7 +196,6 @@ L60:
         if (abse <= thresh) {
             goto L80;
         }
-        smin = min(smin, abss);
         d__1 = max(smax, abss);
         smax = max(d__1, abse);
     }
@@ -243,7 +241,7 @@ L90:
         }
         if (tol >= 0.) {
             mu = (d__1 = d__[ll], abs(d__1));
-            sminl = mu;
+            smin = mu;
             i__1 = m - 1;
             for (lll = ll; lll <= i__1; ++lll) {
                 if ((d__1 = e[lll], abs(d__1)) <= tol * mu) {
@@ -251,7 +249,7 @@ L90:
                     goto L60;
                 }
                 mu = (d__2 = d__[lll + 1], abs(d__2)) * (mu / (mu + (d__1 = e[lll], abs(d__1))));
-                sminl = min(sminl, mu);
+                smin = min(smin, mu);
             }
         }
     } else {
@@ -262,7 +260,7 @@ L90:
         }
         if (tol >= 0.) {
             mu = (d__1 = d__[m], abs(d__1));
-            sminl = mu;
+            smin = mu;
             i__1 = ll;
             for (lll = m - 1; lll >= i__1; --lll) {
                 if ((d__1 = e[lll], abs(d__1)) <= tol * mu) {
@@ -270,14 +268,14 @@ L90:
                     goto L60;
                 }
                 mu = (d__2 = d__[lll], abs(d__2)) * (mu / (mu + (d__1 = e[lll], abs(d__1))));
-                sminl = min(sminl, mu);
+                smin = min(smin, mu);
             }
         }
     }
     oldll = ll;
     oldm = m;
     d__1 = eps, d__2 = tol * .01;
-    if (tol >= 0. && *n * tol * (sminl / smax) <= max(d__1, d__2)) {
+    if (tol >= 0. && *n * tol * (smin / smax) <= max(d__1, d__2)) {
         shift = 0.;
     } else {
         if (idir == 1) {
diff --git a/lib/linalg/dgecon.cpp b/lib/linalg/dgecon.cpp
index 01604f5f5d..71097d622d 100644
--- a/lib/linalg/dgecon.cpp
+++ b/lib/linalg/dgecon.cpp
@@ -20,6 +20,7 @@ int dgecon_(char *norm, integer *n, doublereal *a, integer *lda, doublereal *ano
                 integer *);
     extern doublereal dlamch_(char *, ftnlen);
     extern integer idamax_(integer *, doublereal *, integer *);
+    extern logical disnan_(doublereal *);
     extern int xerbla_(char *, integer *, ftnlen);
     doublereal ainvnm;
     extern int dlatrs_(char *, char *, char *, char *, integer *, doublereal *, integer *,
@@ -27,12 +28,13 @@ int dgecon_(char *norm, integer *n, doublereal *a, integer *lda, doublereal *ano
                        ftnlen);
     logical onenrm;
     char normin[1];
-    doublereal smlnum;
+    doublereal smlnum, hugeval;
     a_dim1 = *lda;
     a_offset = 1 + a_dim1;
     a -= a_offset;
     --work;
     --iwork;
+    hugeval = dlamch_((char *)"Overflow", (ftnlen)8);
     *info = 0;
     onenrm = *(unsigned char *)norm == '1' || lsame_(norm, (char *)"O", (ftnlen)1, (ftnlen)1);
     if (!onenrm && !lsame_(norm, (char *)"I", (ftnlen)1, (ftnlen)1)) {
@@ -55,6 +57,13 @@ int dgecon_(char *norm, integer *n, doublereal *a, integer *lda, doublereal *ano
         return 0;
     } else if (*anorm == 0.) {
         return 0;
+    } else if (disnan_(anorm)) {
+        *rcond = *anorm;
+        *info = -5;
+        return 0;
+    } else if (*anorm > hugeval) {
+        *info = -5;
+        return 0;
     }
     smlnum = dlamch_((char *)"Safe minimum", (ftnlen)12);
     ainvnm = 0.;
@@ -92,6 +101,12 @@ L10:
     }
     if (ainvnm != 0.) {
         *rcond = 1. / ainvnm / *anorm;
+    } else {
+        *info = 1;
+        return 0;
+    }
+    if (disnan_(rcond) || *rcond > hugeval) {
+        *info = 1;
     }
 L20:
     return 0;
diff --git a/lib/linalg/dgelsd.cpp b/lib/linalg/dgelsd.cpp
index 479d95dd61..7f22a11b08 100644
--- a/lib/linalg/dgelsd.cpp
+++ b/lib/linalg/dgelsd.cpp
@@ -19,9 +19,8 @@ int dgelsd_(integer *m, integer *n, integer *nrhs, doublereal *a, integer *lda,
     integer itau, nlvl, iascl, ibscl;
     doublereal sfmin;
     integer minmn, maxmn, itaup, itauq, mnthr, nwork;
-    extern int dlabad_(doublereal *, doublereal *),
-        dgebrd_(integer *, integer *, doublereal *, integer *, doublereal *, doublereal *,
-                doublereal *, doublereal *, doublereal *, integer *, integer *);
+    extern int dgebrd_(integer *, integer *, doublereal *, integer *, doublereal *, doublereal *,
+                       doublereal *, doublereal *, doublereal *, integer *, integer *);
     extern doublereal dlamch_(char *, ftnlen),
         dlange_(char *, integer *, integer *, doublereal *, integer *, doublereal *, ftnlen);
     extern int dgelqf_(integer *, integer *, doublereal *, integer *, doublereal *, doublereal *,
@@ -189,7 +188,6 @@ int dgelsd_(integer *m, integer *n, integer *nrhs, doublereal *a, integer *lda,
     sfmin = dlamch_((char *)"S", (ftnlen)1);
     smlnum = sfmin / eps;
     bignum = 1. / smlnum;
-    dlabad_(&smlnum, &bignum);
     anrm = dlange_((char *)"M", m, n, &a[a_offset], lda, &work[1], (ftnlen)1);
     iascl = 0;
     if (anrm > 0. && anrm < smlnum) {
diff --git a/lib/linalg/dgelss.cpp b/lib/linalg/dgelss.cpp
index e10906f4e9..a11ea7e297 100644
--- a/lib/linalg/dgelss.cpp
+++ b/lib/linalg/dgelss.cpp
@@ -30,9 +30,8 @@ int dgelss_(integer *m, integer *n, integer *nrhs, doublereal *a, integer *lda,
     integer minmn;
     extern int dcopy_(integer *, doublereal *, integer *, doublereal *, integer *);
     integer maxmn, itaup, itauq, mnthr, iwork;
-    extern int dlabad_(doublereal *, doublereal *),
-        dgebrd_(integer *, integer *, doublereal *, integer *, doublereal *, doublereal *,
-                doublereal *, doublereal *, doublereal *, integer *, integer *);
+    extern int dgebrd_(integer *, integer *, doublereal *, integer *, doublereal *, doublereal *,
+                       doublereal *, doublereal *, doublereal *, integer *, integer *);
     extern doublereal dlamch_(char *, ftnlen),
         dlange_(char *, integer *, integer *, doublereal *, integer *, doublereal *, ftnlen);
     integer bdspac;
@@ -208,7 +207,6 @@ int dgelss_(integer *m, integer *n, integer *nrhs, doublereal *a, integer *lda,
     sfmin = dlamch_((char *)"S", (ftnlen)1);
     smlnum = sfmin / eps;
     bignum = 1. / smlnum;
-    dlabad_(&smlnum, &bignum);
     anrm = dlange_((char *)"M", m, n, &a[a_offset], lda, &work[1], (ftnlen)1);
     iascl = 0;
     if (anrm > 0. && anrm < smlnum) {
@@ -300,7 +298,7 @@ int dgelss_(integer *m, integer *n, integer *nrhs, doublereal *a, integer *lda,
                        &c_b46, &work[1], n, (ftnlen)1, (ftnlen)1);
                 dlacpy_((char *)"G", n, &bl, &work[1], n, &b[i__ * b_dim1 + 1], ldb, (ftnlen)1);
             }
-        } else {
+        } else if (*nrhs == 1) {
             dgemv_((char *)"T", n, n, &c_b79, &a[a_offset], lda, &b[b_offset], &c__1, &c_b46, &work[1],
                    &c__1, (ftnlen)1);
             dcopy_(n, &work[1], &c__1, &b[b_offset], &c__1);
@@ -376,7 +374,7 @@ int dgelss_(integer *m, integer *n, integer *nrhs, doublereal *a, integer *lda,
                            ldb, &c_b46, &work[iwork], m, (ftnlen)1, (ftnlen)1);
                     dlacpy_((char *)"G", m, &bl, &work[iwork], m, &b[i__ * b_dim1 + 1], ldb, (ftnlen)1);
                 }
-            } else {
+            } else if (*nrhs == 1) {
                 dgemv_((char *)"T", m, m, &c_b79, &work[il], &ldwork, &b[b_dim1 + 1], &c__1, &c_b46,
                        &work[iwork], &c__1, (ftnlen)1);
                 dcopy_(m, &work[iwork], &c__1, &b[b_dim1 + 1], &c__1);
@@ -438,7 +436,7 @@ int dgelss_(integer *m, integer *n, integer *nrhs, doublereal *a, integer *lda,
                            ldb, &c_b46, &work[1], n, (ftnlen)1, (ftnlen)1);
                     dlacpy_((char *)"F", n, &bl, &work[1], n, &b[i__ * b_dim1 + 1], ldb, (ftnlen)1);
                 }
-            } else {
+            } else if (*nrhs == 1) {
                 dgemv_((char *)"T", m, n, &c_b79, &a[a_offset], lda, &b[b_offset], &c__1, &c_b46, &work[1],
                        &c__1, (ftnlen)1);
                 dcopy_(n, &work[1], &c__1, &b[b_offset], &c__1);
diff --git a/lib/linalg/dlabad.cpp b/lib/linalg/dlabad.cpp
index 790163be3c..760cfd8f20 100644
--- a/lib/linalg/dlabad.cpp
+++ b/lib/linalg/dlabad.cpp
@@ -4,11 +4,6 @@ extern "C" {
 #include "lmp_f2c.h"
 int dlabad_(doublereal *small, doublereal *large)
 {
-    double d_lmp_lg10(doublereal *), sqrt(doublereal);
-    if (d_lmp_lg10(large) > 2e3) {
-        *small = sqrt(*small);
-        *large = sqrt(*large);
-    }
     return 0;
 }
 #ifdef __cplusplus
diff --git a/lib/linalg/dlaed2.cpp b/lib/linalg/dlaed2.cpp
index 4f0461bab9..fddf272a82 100644
--- a/lib/linalg/dlaed2.cpp
+++ b/lib/linalg/dlaed2.cpp
@@ -5,7 +5,7 @@ extern "C" {
 static doublereal c_b3 = -1.;
 static integer c__1 = 1;
 int dlaed2_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q, integer *ldq,
-            integer *indxq, doublereal *rho, doublereal *z__, doublereal *dlamda, doublereal *w,
+            integer *indxq, doublereal *rho, doublereal *z__, doublereal *dlambda, doublereal *w,
             doublereal *q2, integer *indx, integer *indxc, integer *indxp, integer *coltyp,
             integer *info)
 {
@@ -35,7 +35,7 @@ int dlaed2_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q,
     q -= q_offset;
     --indxq;
     --z__;
-    --dlamda;
+    --dlambda;
     --w;
     --q2;
     --indx;
@@ -75,9 +75,9 @@ int dlaed2_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q,
     }
     i__1 = *n;
     for (i__ = 1; i__ <= i__1; ++i__) {
-        dlamda[i__] = d__[indxq[i__]];
+        dlambda[i__] = d__[indxq[i__]];
     }
-    dlamrg_(n1, &n2, &dlamda[1], &c__1, &c__1, &indxc[1]);
+    dlamrg_(n1, &n2, &dlambda[1], &c__1, &c__1, &indxc[1]);
     i__1 = *n;
     for (i__ = 1; i__ <= i__1; ++i__) {
         indx[i__] = indxq[indxc[i__]];
@@ -94,11 +94,11 @@ int dlaed2_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q,
         for (j = 1; j <= i__1; ++j) {
             i__ = indx[j];
             dcopy_(n, &q[i__ * q_dim1 + 1], &c__1, &q2[iq2], &c__1);
-            dlamda[j] = d__[i__];
+            dlambda[j] = d__[i__];
             iq2 += *n;
         }
         dlacpy_((char *)"A", n, n, &q2[1], n, &q[q_offset], ldq, (ftnlen)1);
-        dcopy_(n, &dlamda[1], &c__1, &d__[1], &c__1);
+        dcopy_(n, &dlambda[1], &c__1, &d__[1], &c__1);
         goto L190;
     }
     i__1 = *n1;
@@ -176,7 +176,7 @@ L80:
             pj = nj;
         } else {
             ++(*k);
-            dlamda[*k] = d__[pj];
+            dlambda[*k] = d__[pj];
             w[*k] = z__[pj];
             indxp[*k] = pj;
             pj = nj;
@@ -185,7 +185,7 @@ L80:
     goto L80;
 L100:
     ++(*k);
-    dlamda[*k] = d__[pj];
+    dlambda[*k] = d__[pj];
     w[*k] = z__[pj];
     indxp[*k] = pj;
     for (j = 1; j <= 4; ++j) {
diff --git a/lib/linalg/dlaed3.cpp b/lib/linalg/dlaed3.cpp
index 926b0ecd7a..1c08f7d3de 100644
--- a/lib/linalg/dlaed3.cpp
+++ b/lib/linalg/dlaed3.cpp
@@ -3,10 +3,10 @@ extern "C" {
 #endif
 #include "lmp_f2c.h"
 static integer c__1 = 1;
-static doublereal c_b22 = 1.;
-static doublereal c_b23 = 0.;
+static doublereal c_b21 = 1.;
+static doublereal c_b22 = 0.;
 int dlaed3_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q, integer *ldq,
-            doublereal *rho, doublereal *dlamda, doublereal *q2, integer *indx, integer *ctot,
+            doublereal *rho, doublereal *dlambda, doublereal *q2, integer *indx, integer *ctot,
             doublereal *w, doublereal *s, integer *info)
 {
     integer q_dim1, q_offset, i__1, i__2;
@@ -20,10 +20,9 @@ int dlaed3_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q,
                       ftnlen, ftnlen),
         dcopy_(integer *, doublereal *, integer *, doublereal *, integer *),
         dlaed4_(integer *, integer *, doublereal *, doublereal *, doublereal *, doublereal *,
-                doublereal *, integer *);
-    extern doublereal dlamc3_(doublereal *, doublereal *);
-    extern int dlacpy_(char *, integer *, integer *, doublereal *, integer *, doublereal *,
-                       integer *, ftnlen),
+                doublereal *, integer *),
+        dlacpy_(char *, integer *, integer *, doublereal *, integer *, doublereal *, integer *,
+                ftnlen),
         dlaset_(char *, integer *, integer *, doublereal *, doublereal *, doublereal *, integer *,
                 ftnlen),
         xerbla_(char *, integer *, ftnlen);
@@ -31,7 +30,7 @@ int dlaed3_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q,
     q_dim1 = *ldq;
     q_offset = 1 + q_dim1;
     q -= q_offset;
-    --dlamda;
+    --dlambda;
     --q2;
     --indx;
     --ctot;
@@ -54,12 +53,8 @@ int dlaed3_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q,
         return 0;
     }
     i__1 = *k;
-    for (i__ = 1; i__ <= i__1; ++i__) {
-        dlamda[i__] = dlamc3_(&dlamda[i__], &dlamda[i__]) - dlamda[i__];
-    }
-    i__1 = *k;
     for (j = 1; j <= i__1; ++j) {
-        dlaed4_(k, &j, &dlamda[1], &w[1], &q[j * q_dim1 + 1], rho, &d__[j], info);
+        dlaed4_(k, &j, &dlambda[1], &w[1], &q[j * q_dim1 + 1], rho, &d__[j], info);
         if (*info != 0) {
             goto L120;
         }
@@ -86,11 +81,11 @@ int dlaed3_(integer *k, integer *n, integer *n1, doublereal *d__, doublereal *q,
     for (j = 1; j <= i__1; ++j) {
         i__2 = j - 1;
         for (i__ = 1; i__ <= i__2; ++i__) {
-            w[i__] *= q[i__ + j * q_dim1] / (dlamda[i__] - dlamda[j]);
+            w[i__] *= q[i__ + j * q_dim1] / (dlambda[i__] - dlambda[j]);
         }
         i__2 = *k;
         for (i__ = j + 1; i__ <= i__2; ++i__) {
-            w[i__] *= q[i__ + j * q_dim1] / (dlamda[i__] - dlamda[j]);
+            w[i__] *= q[i__ + j * q_dim1] / (dlambda[i__] - dlambda[j]);
         }
     }
     i__1 = *k;
@@ -118,17 +113,17 @@ L110:
     dlacpy_((char *)"A", &n23, k, &q[ctot[1] + 1 + q_dim1], ldq, &s[1], &n23, (ftnlen)1);
     iq2 = *n1 * n12 + 1;
     if (n23 != 0) {
-        dgemm_((char *)"N", (char *)"N", &n2, k, &n23, &c_b22, &q2[iq2], &n2, &s[1], &n23, &c_b23,
+        dgemm_((char *)"N", (char *)"N", &n2, k, &n23, &c_b21, &q2[iq2], &n2, &s[1], &n23, &c_b22,
                &q[*n1 + 1 + q_dim1], ldq, (ftnlen)1, (ftnlen)1);
     } else {
-        dlaset_((char *)"A", &n2, k, &c_b23, &c_b23, &q[*n1 + 1 + q_dim1], ldq, (ftnlen)1);
+        dlaset_((char *)"A", &n2, k, &c_b22, &c_b22, &q[*n1 + 1 + q_dim1], ldq, (ftnlen)1);
     }
     dlacpy_((char *)"A", &n12, k, &q[q_offset], ldq, &s[1], &n12, (ftnlen)1);
     if (n12 != 0) {
-        dgemm_((char *)"N", (char *)"N", n1, k, &n12, &c_b22, &q2[1], n1, &s[1], &n12, &c_b23, &q[q_offset], ldq,
+        dgemm_((char *)"N", (char *)"N", n1, k, &n12, &c_b21, &q2[1], n1, &s[1], &n12, &c_b22, &q[q_offset], ldq,
                (ftnlen)1, (ftnlen)1);
     } else {
-        dlaset_((char *)"A", n1, k, &c_b23, &c_b23, &q[q_dim1 + 1], ldq, (ftnlen)1);
+        dlaset_((char *)"A", n1, k, &c_b22, &c_b22, &q[q_dim1 + 1], ldq, (ftnlen)1);
     }
 L120:
     return 0;
diff --git a/lib/linalg/dlaed8.cpp b/lib/linalg/dlaed8.cpp
index 46580ce44f..088d4a3e9d 100644
--- a/lib/linalg/dlaed8.cpp
+++ b/lib/linalg/dlaed8.cpp
@@ -6,7 +6,7 @@ static doublereal c_b3 = -1.;
 static integer c__1 = 1;
 int dlaed8_(integer *icompq, integer *k, integer *n, integer *qsiz, doublereal *d__, doublereal *q,
             integer *ldq, integer *indxq, doublereal *rho, integer *cutpnt, doublereal *z__,
-            doublereal *dlamda, doublereal *q2, integer *ldq2, doublereal *w, integer *perm,
+            doublereal *dlambda, doublereal *q2, integer *ldq2, doublereal *w, integer *perm,
             integer *givptr, integer *givcol, doublereal *givnum, integer *indxp, integer *indx,
             integer *info)
 {
@@ -35,7 +35,7 @@ int dlaed8_(integer *icompq, integer *k, integer *n, integer *qsiz, doublereal *
     q -= q_offset;
     --indxq;
     --z__;
-    --dlamda;
+    --dlambda;
     q2_dim1 = *ldq2;
     q2_offset = 1 + q2_dim1;
     q2 -= q2_offset;
@@ -87,15 +87,15 @@ int dlaed8_(integer *icompq, integer *k, integer *n, integer *qsiz, doublereal *
     }
     i__1 = *n;
     for (i__ = 1; i__ <= i__1; ++i__) {
-        dlamda[i__] = d__[indxq[i__]];
+        dlambda[i__] = d__[indxq[i__]];
         w[i__] = z__[indxq[i__]];
     }
     i__ = 1;
     j = *cutpnt + 1;
-    dlamrg_(&n1, &n2, &dlamda[1], &c__1, &c__1, &indx[1]);
+    dlamrg_(&n1, &n2, &dlambda[1], &c__1, &c__1, &indx[1]);
     i__1 = *n;
     for (i__ = 1; i__ <= i__1; ++i__) {
-        d__[i__] = dlamda[indx[i__]];
+        d__[i__] = dlambda[indx[i__]];
         z__[i__] = w[indx[i__]];
     }
     imax = idamax_(n, &z__[1], &c__1);
@@ -183,7 +183,7 @@ L80:
         } else {
             ++(*k);
             w[*k] = z__[jlam];
-            dlamda[*k] = d__[jlam];
+            dlambda[*k] = d__[jlam];
             indxp[*k] = jlam;
             jlam = j;
         }
@@ -192,21 +192,21 @@ L80:
 L100:
     ++(*k);
     w[*k] = z__[jlam];
-    dlamda[*k] = d__[jlam];
+    dlambda[*k] = d__[jlam];
     indxp[*k] = jlam;
 L110:
     if (*icompq == 0) {
         i__1 = *n;
         for (j = 1; j <= i__1; ++j) {
             jp = indxp[j];
-            dlamda[j] = d__[jp];
+            dlambda[j] = d__[jp];
             perm[j] = indxq[indx[jp]];
         }
     } else {
         i__1 = *n;
         for (j = 1; j <= i__1; ++j) {
             jp = indxp[j];
-            dlamda[j] = d__[jp];
+            dlambda[j] = d__[jp];
             perm[j] = indxq[indx[jp]];
             dcopy_(qsiz, &q[perm[j] * q_dim1 + 1], &c__1, &q2[j * q2_dim1 + 1], &c__1);
         }
@@ -214,10 +214,10 @@ L110:
     if (*k < *n) {
         if (*icompq == 0) {
             i__1 = *n - *k;
-            dcopy_(&i__1, &dlamda[*k + 1], &c__1, &d__[*k + 1], &c__1);
+            dcopy_(&i__1, &dlambda[*k + 1], &c__1, &d__[*k + 1], &c__1);
         } else {
             i__1 = *n - *k;
-            dcopy_(&i__1, &dlamda[*k + 1], &c__1, &d__[*k + 1], &c__1);
+            dcopy_(&i__1, &dlambda[*k + 1], &c__1, &d__[*k + 1], &c__1);
             i__1 = *n - *k;
             dlacpy_((char *)"A", qsiz, &i__1, &q2[(*k + 1) * q2_dim1 + 1], ldq2, &q[(*k + 1) * q_dim1 + 1],
                     ldq, (ftnlen)1);
diff --git a/lib/linalg/dlaed9.cpp b/lib/linalg/dlaed9.cpp
index 2ca15ee0d7..76aa3bc747 100644
--- a/lib/linalg/dlaed9.cpp
+++ b/lib/linalg/dlaed9.cpp
@@ -4,7 +4,7 @@ extern "C" {
 #include "lmp_f2c.h"
 static integer c__1 = 1;
 int dlaed9_(integer *k, integer *kstart, integer *kstop, integer *n, doublereal *d__, doublereal *q,
-            integer *ldq, doublereal *rho, doublereal *dlamda, doublereal *w, doublereal *s,
+            integer *ldq, doublereal *rho, doublereal *dlambda, doublereal *w, doublereal *s,
             integer *lds, integer *info)
 {
     integer q_dim1, q_offset, s_dim1, s_offset, i__1, i__2;
@@ -15,14 +15,13 @@ int dlaed9_(integer *k, integer *kstart, integer *kstop, integer *n, doublereal
     extern doublereal dnrm2_(integer *, doublereal *, integer *);
     extern int dcopy_(integer *, doublereal *, integer *, doublereal *, integer *),
         dlaed4_(integer *, integer *, doublereal *, doublereal *, doublereal *, doublereal *,
-                doublereal *, integer *);
-    extern doublereal dlamc3_(doublereal *, doublereal *);
-    extern int xerbla_(char *, integer *, ftnlen);
+                doublereal *, integer *),
+        xerbla_(char *, integer *, ftnlen);
     --d__;
     q_dim1 = *ldq;
     q_offset = 1 + q_dim1;
     q -= q_offset;
-    --dlamda;
+    --dlambda;
     --w;
     s_dim1 = *lds;
     s_offset = 1 + s_dim1;
@@ -49,13 +48,9 @@ int dlaed9_(integer *k, integer *kstart, integer *kstop, integer *n, doublereal
     if (*k == 0) {
         return 0;
     }
-    i__1 = *n;
-    for (i__ = 1; i__ <= i__1; ++i__) {
-        dlamda[i__] = dlamc3_(&dlamda[i__], &dlamda[i__]) - dlamda[i__];
-    }
     i__1 = *kstop;
     for (j = *kstart; j <= i__1; ++j) {
-        dlaed4_(k, &j, &dlamda[1], &w[1], &q[j * q_dim1 + 1], rho, &d__[j], info);
+        dlaed4_(k, &j, &dlambda[1], &w[1], &q[j * q_dim1 + 1], rho, &d__[j], info);
         if (*info != 0) {
             goto L120;
         }
@@ -77,11 +72,11 @@ int dlaed9_(integer *k, integer *kstart, integer *kstop, integer *n, doublereal
     for (j = 1; j <= i__1; ++j) {
         i__2 = j - 1;
         for (i__ = 1; i__ <= i__2; ++i__) {
-            w[i__] *= q[i__ + j * q_dim1] / (dlamda[i__] - dlamda[j]);
+            w[i__] *= q[i__ + j * q_dim1] / (dlambda[i__] - dlambda[j]);
         }
         i__2 = *k;
         for (i__ = j + 1; i__ <= i__2; ++i__) {
-            w[i__] *= q[i__ + j * q_dim1] / (dlamda[i__] - dlamda[j]);
+            w[i__] *= q[i__ + j * q_dim1] / (dlambda[i__] - dlambda[j]);
         }
     }
     i__1 = *k;
diff --git a/lib/linalg/dlasd8.cpp b/lib/linalg/dlasd8.cpp
index ba92932435..57f353ddd4 100644
--- a/lib/linalg/dlasd8.cpp
+++ b/lib/linalg/dlasd8.cpp
@@ -4,7 +4,7 @@ extern "C" {
 #include "lmp_f2c.h"
 static integer c__1 = 1;
 static integer c__0 = 0;
-static doublereal c_b8 = 1.;
+static doublereal c_b7 = 1.;
 int dlasd8_(integer *icompq, integer *k, doublereal *d__, doublereal *z__, doublereal *vf,
             doublereal *vl, doublereal *difl, doublereal *difr, integer *lddifr, doublereal *dsigma,
             doublereal *work, integer *info)
@@ -62,19 +62,15 @@ int dlasd8_(integer *icompq, integer *k, doublereal *d__, doublereal *z__, doubl
         }
         return 0;
     }
-    i__1 = *k;
-    for (i__ = 1; i__ <= i__1; ++i__) {
-        dsigma[i__] = dlamc3_(&dsigma[i__], &dsigma[i__]) - dsigma[i__];
-    }
     iwk1 = 1;
     iwk2 = iwk1 + *k;
     iwk3 = iwk2 + *k;
     iwk2i = iwk2 - 1;
     iwk3i = iwk3 - 1;
     rho = dnrm2_(k, &z__[1], &c__1);
-    dlascl_((char *)"G", &c__0, &c__0, &rho, &c_b8, k, &c__1, &z__[1], k, info, (ftnlen)1);
+    dlascl_((char *)"G", &c__0, &c__0, &rho, &c_b7, k, &c__1, &z__[1], k, info, (ftnlen)1);
     rho *= rho;
-    dlaset_((char *)"A", k, &c__1, &c_b8, &c_b8, &work[iwk3], k, (ftnlen)1);
+    dlaset_((char *)"A", k, &c__1, &c_b7, &c_b7, &work[iwk3], k, (ftnlen)1);
     i__1 = *k;
     for (j = 1; j <= i__1; ++j) {
         dlasd4_(k, &j, &dsigma[1], &z__[1], &work[iwk1], &rho, &d__[j], &work[iwk2], info);
diff --git a/lib/linalg/dlatrs.cpp b/lib/linalg/dlatrs.cpp
index bd2af669dc..b7af6b32a9 100644
--- a/lib/linalg/dlatrs.cpp
+++ b/lib/linalg/dlatrs.cpp
@@ -16,7 +16,7 @@ int dlatrs_(char *uplo, char *trans, char *diag, char *normin, integer *n, doubl
     extern doublereal ddot_(integer *, doublereal *, integer *, doublereal *, integer *);
     doublereal xbnd;
     integer imax;
-    doublereal tmax, tjjs, xmax, grow, sumj;
+    doublereal tmax, tjjs, xmax, grow, sumj, work[1];
     extern int dscal_(integer *, doublereal *, doublereal *, integer *);
     extern logical lsame_(char *, char *, ftnlen, ftnlen);
     doublereal tscal, uscal;
@@ -100,7 +100,7 @@ int dlatrs_(char *uplo, char *trans, char *diag, char *normin, integer *n, doubl
                 i__1 = *n;
                 for (j = 2; j <= i__1; ++j) {
                     i__2 = j - 1;
-                    d__1 = dlange_((char *)"M", &i__2, &c__1, &a[j * a_dim1 + 1], &c__1, &sumj, (ftnlen)1);
+                    d__1 = dlange_((char *)"M", &i__2, &c__1, &a[j * a_dim1 + 1], &c__1, work, (ftnlen)1);
                     tmax = max(d__1, tmax);
                 }
             } else {
@@ -108,7 +108,7 @@ int dlatrs_(char *uplo, char *trans, char *diag, char *normin, integer *n, doubl
                 for (j = 1; j <= i__1; ++j) {
                     i__2 = *n - j;
                     d__1 =
-                        dlange_((char *)"M", &i__2, &c__1, &a[j + 1 + j * a_dim1], &c__1, &sumj, (ftnlen)1);
+                        dlange_((char *)"M", &i__2, &c__1, &a[j + 1 + j * a_dim1], &c__1, work, (ftnlen)1);
                     tmax = max(d__1, tmax);
                 }
             }
diff --git a/lib/linalg/drscl.cpp b/lib/linalg/drscl.cpp
index 90e278a709..bd35c4b276 100644
--- a/lib/linalg/drscl.cpp
+++ b/lib/linalg/drscl.cpp
@@ -7,8 +7,7 @@ int drscl_(integer *n, doublereal *sa, doublereal *sx, integer *incx)
     doublereal mul, cden;
     logical done;
     doublereal cnum, cden1, cnum1;
-    extern int dscal_(integer *, doublereal *, doublereal *, integer *),
-        dlabad_(doublereal *, doublereal *);
+    extern int dscal_(integer *, doublereal *, doublereal *, integer *);
     extern doublereal dlamch_(char *, ftnlen);
     doublereal bignum, smlnum;
     --sx;
@@ -17,7 +16,6 @@ int drscl_(integer *n, doublereal *sa, doublereal *sx, integer *incx)
     }
     smlnum = dlamch_((char *)"S", (ftnlen)1);
     bignum = 1. / smlnum;
-    dlabad_(&smlnum, &bignum);
     cden = *sa;
     cnum = 1.;
 L10:
diff --git a/lib/linalg/ilaenv.cpp b/lib/linalg/ilaenv.cpp
index 1cc1c571f1..4f5a65a6ea 100644
--- a/lib/linalg/ilaenv.cpp
+++ b/lib/linalg/ilaenv.cpp
@@ -3,8 +3,8 @@ extern "C" {
 #endif
 #include "lmp_f2c.h"
 static integer c__1 = 1;
-static real c_b176 = (float)0.;
-static real c_b177 = (float)1.;
+static real c_b179 = (float)0.;
+static real c_b180 = (float)1.;
 static integer c__0 = 0;
 integer ilaenv_(integer *ispec, char *name__, char *opts, integer *n1, integer *n2, integer *n3,
                 integer *n4, ftnlen name_len, ftnlen opts_len)
@@ -201,6 +201,12 @@ L50:
             } else {
                 nb = 64;
             }
+        } else if (s_lmp_cmp(subnam + 3, (char *)"QP3RK", (ftnlen)4, (ftnlen)5) == 0) {
+            if (sname) {
+                nb = 32;
+            } else {
+                nb = 32;
+            }
         }
     } else if (s_lmp_cmp(c2, (char *)"PO", (ftnlen)2, (ftnlen)2) == 0) {
         if (s_lmp_cmp(c3, (char *)"TRF", (ftnlen)3, (ftnlen)3) == 0) {
@@ -402,6 +408,12 @@ L60:
             } else {
                 nbmin = 2;
             }
+        } else if (s_lmp_cmp(subnam + 3, (char *)"QP3RK", (ftnlen)4, (ftnlen)5) == 0) {
+            if (sname) {
+                nbmin = 2;
+            } else {
+                nbmin = 2;
+            }
         }
     } else if (s_lmp_cmp(c2, (char *)"SY", (ftnlen)2, (ftnlen)2) == 0) {
         if (s_lmp_cmp(c3, (char *)"TRF", (ftnlen)3, (ftnlen)3) == 0) {
@@ -493,6 +505,12 @@ L70:
             } else {
                 nx = 128;
             }
+        } else if (s_lmp_cmp(subnam + 3, (char *)"QP3RK", (ftnlen)4, (ftnlen)5) == 0) {
+            if (sname) {
+                nx = 128;
+            } else {
+                nx = 128;
+            }
         }
     } else if (s_lmp_cmp(c2, (char *)"SY", (ftnlen)2, (ftnlen)2) == 0) {
         if (sname && s_lmp_cmp(c3, (char *)"TRD", (ftnlen)3, (ftnlen)3) == 0) {
@@ -555,13 +573,13 @@ L130:
 L140:
     ret_val = 1;
     if (ret_val == 1) {
-        ret_val = ieeeck_(&c__1, &c_b176, &c_b177);
+        ret_val = ieeeck_(&c__1, &c_b179, &c_b180);
     }
     return ret_val;
 L150:
     ret_val = 1;
     if (ret_val == 1) {
-        ret_val = ieeeck_(&c__0, &c_b176, &c_b177);
+        ret_val = ieeeck_(&c__0, &c_b179, &c_b180);
     }
     return ret_val;
 L160:
diff --git a/lib/linalg/zlaed8.cpp b/lib/linalg/zlaed8.cpp
index d29b380587..0c06006ecd 100644
--- a/lib/linalg/zlaed8.cpp
+++ b/lib/linalg/zlaed8.cpp
@@ -5,7 +5,7 @@ extern "C" {
 static doublereal c_b3 = -1.;
 static integer c__1 = 1;
 int zlaed8_(integer *k, integer *n, integer *qsiz, doublecomplex *q, integer *ldq, doublereal *d__,
-            doublereal *rho, integer *cutpnt, doublereal *z__, doublereal *dlamda,
+            doublereal *rho, integer *cutpnt, doublereal *z__, doublereal *dlambda,
             doublecomplex *q2, integer *ldq2, doublereal *w, integer *indxp, integer *indx,
             integer *indxq, integer *perm, integer *givptr, integer *givcol, doublereal *givnum,
             integer *info)
@@ -35,7 +35,7 @@ int zlaed8_(integer *k, integer *n, integer *qsiz, doublecomplex *q, integer *ld
     q -= q_offset;
     --d__;
     --z__;
-    --dlamda;
+    --dlambda;
     q2_dim1 = *ldq2;
     q2_offset = 1 + q2_dim1;
     q2 -= q2_offset;
@@ -86,15 +86,15 @@ int zlaed8_(integer *k, integer *n, integer *qsiz, doublecomplex *q, integer *ld
     }
     i__1 = *n;
     for (i__ = 1; i__ <= i__1; ++i__) {
-        dlamda[i__] = d__[indxq[i__]];
+        dlambda[i__] = d__[indxq[i__]];
         w[i__] = z__[indxq[i__]];
     }
     i__ = 1;
     j = *cutpnt + 1;
-    dlamrg_(&n1, &n2, &dlamda[1], &c__1, &c__1, &indx[1]);
+    dlamrg_(&n1, &n2, &dlambda[1], &c__1, &c__1, &indx[1]);
     i__1 = *n;
     for (i__ = 1; i__ <= i__1; ++i__) {
-        d__[i__] = dlamda[indx[i__]];
+        d__[i__] = dlambda[indx[i__]];
         z__[i__] = w[indx[i__]];
     }
     imax = idamax_(n, &z__[1], &c__1);
@@ -173,7 +173,7 @@ L70:
         } else {
             ++(*k);
             w[*k] = z__[jlam];
-            dlamda[*k] = d__[jlam];
+            dlambda[*k] = d__[jlam];
             indxp[*k] = jlam;
             jlam = j;
         }
@@ -182,19 +182,19 @@ L70:
 L90:
     ++(*k);
     w[*k] = z__[jlam];
-    dlamda[*k] = d__[jlam];
+    dlambda[*k] = d__[jlam];
     indxp[*k] = jlam;
 L100:
     i__1 = *n;
     for (j = 1; j <= i__1; ++j) {
         jp = indxp[j];
-        dlamda[j] = d__[jp];
+        dlambda[j] = d__[jp];
         perm[j] = indxq[indx[jp]];
         zcopy_(qsiz, &q[perm[j] * q_dim1 + 1], &c__1, &q2[j * q2_dim1 + 1], &c__1);
     }
     if (*k < *n) {
         i__1 = *n - *k;
-        dcopy_(&i__1, &dlamda[*k + 1], &c__1, &d__[*k + 1], &c__1);
+        dcopy_(&i__1, &dlambda[*k + 1], &c__1, &d__[*k + 1], &c__1);
         i__1 = *n - *k;
         zlacpy_((char *)"A", qsiz, &i__1, &q2[(*k + 1) * q2_dim1 + 1], ldq2, &q[(*k + 1) * q_dim1 + 1], ldq,
                 (ftnlen)1);