diff --git a/src/atmesc.c b/src/atmesc.c
index ee34af7d3..62d8a66f1 100644
--- a/src/atmesc.c
+++ b/src/atmesc.c
@@ -1692,6 +1692,22 @@ void fnForceBehaviorAtmEsc(BODY *body, MODULE *module, EVOLVE *evolve, IO *io,
   }
 }
 
+void AuxPropsLehmer17(BODY *body,int iBody) {
+  if (body[iBody].bAutoThermTemp) {
+      body[iBody].dThermTemp = fdThermalTemp(body, iBody);
+    }
+    body[iBody].dGravAccel = BIGG *
+                             (body[iBody].dMass - body[iBody].dEnvelopeMass) /
+                             (body[iBody].dRadSolid * body[iBody].dRadSolid);
+    body[iBody].dScaleHeight = body[iBody].dAtmGasConst *
+                               body[iBody].dThermTemp / body[iBody].dGravAccel;
+    body[iBody].dPresSurf =
+          fdLehmerPres(body[iBody].dEnvelopeMass, body[iBody].dGravAccel,
+                       body[iBody].dRadSolid);
+    body[iBody].dRadXUV = fdLehmerRadius(body, iBody);
+    body[iBody].dRadius = body[iBody].dRadXUV / body[iBody].dXFrac;
+}
+
 
 /**
 Initializes several helper variables and properties used in the integration.
@@ -1707,41 +1723,13 @@ Initializes several helper variables and properties used in the integration.
 void fnPropsAuxAtmEsc(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update,
                       int iBody) {
 
-  /*
-  #ifdef DEBUG
-    if (body[iBody].dMass < 0) {
-      fprintf(stderr,"ERROR: %s's mass is %.5e at %.5e
-  years.\n",body[iBody].cName, body[iBody].dMass,evolve->dTime/YEARSEC);
-      exit(EXIT_INT);
-    }
-  #endif
-  */
-
   if (body[iBody].iPlanetRadiusModel == ATMESC_LEHMER17) {
-    if (body[iBody].bAutoThermTemp) {
-      body[iBody].dThermTemp = fdThermalTemp(body, iBody);
-    }
-    body[iBody].dGravAccel = BIGG *
-                             (body[iBody].dMass - body[iBody].dEnvelopeMass) /
-                             (body[iBody].dRadSolid * body[iBody].dRadSolid);
-    body[iBody].dScaleHeight = body[iBody].dAtmGasConst *
-                               body[iBody].dThermTemp / body[iBody].dGravAccel;
-    body[iBody].dPresSurf =
-          fdLehmerPres(body[iBody].dEnvelopeMass, body[iBody].dGravAccel,
-                       body[iBody].dRadSolid);
-    body[iBody].dRadXUV = fdLehmerRadius(body, iBody);
-    body[iBody].dRadius = body[iBody].dRadXUV / body[iBody].dXFrac;
+    AuxPropsLehmer17(body,iBody);
   }
 
   // Compute various radii of interest
   body[iBody].dBondiRadius = fdBondiRadius(body, iBody);
   body[iBody].dRocheRadius = fdRocheRadius(body, iBody);
-
-  // Ktide (due to body zero only). WARNING: not suited for binary...
-  /* xi = (pow(body[iBody].dMass / (3. * body[0].dMass), (1. / 3)) *
-               body[iBody].dSemi) /
-              (body[iBody].dRadius * body[iBody].dXFrac);
-  */
   body[iBody].dAtmEscXi = fdAtmEscXi(body, iBody);
   body[iBody].dKTide    = fdKTide(body, io, iBody);
 
diff --git a/src/body.c b/src/body.c
index a3f0dd35e..736b2bc0f 100644
--- a/src/body.c
+++ b/src/body.c
@@ -665,16 +665,16 @@ double fdImK2Total(BODY *body, int iBody) {
   @return Upper mantle k2 Love number
 */
 double fdK2Man(BODY *body, int iBody) {
-  //  return 1.5/(1+9.5*body[iBody].dShmodUMan/(STIFFNESS));
-  return 1.5 / (1 + 9.5 * body[iBody].dShmodUMan / body[iBody].dStiffness);
+  double dK2Man = 1.5 / (1 + 9.5 * body[iBody].dShmodUMan / body[iBody].dStiffness);
+  return dK2Man;
 }
 
 double fdTidalQMan(BODY *body, int iBody) {
   double ShmodUManArr;
   ShmodUManArr = body[iBody].dShmodUMan * body[iBody].dMeltfactorUMan;
-  // return body[iBody].dDynamViscos*body[iBody].dMeanMotion/ShmodUManArr;
-  return body[iBody].dDynamViscos * body[iBody].dMeanMotion /
+  double dTidalQMan = body[iBody].dDynamViscos * body[iBody].dMeanMotion /
          body[iBody].dShmodUMan;
+  return dTidalQMan;
 }
 
 double fdImK2ManThermint(BODY *body, int iBody) {
diff --git a/src/distorb.c b/src/distorb.c
index e6497ff58..f9b000a09 100644
--- a/src/distorb.c
+++ b/src/distorb.c
@@ -1645,7 +1645,7 @@ int fniHaltHillStab(BODY *body, EVOLVE *evolve, HALT *halt, IO *io,
           gamma2 = sqrt(1 - (body[iBody].dHecc * body[iBody].dHecc +
                              body[iBody].dKecc * body[iBody].dKecc));
           delta  = sqrt(body[iBody].dSemi / body[jBody].dSemi);
-        } else if (body[jBody].dSemi > body[iBody].dSemi) {
+        } else {
           // jBody is the outer planet
           mu2    = body[jBody].dMass / body[0].dMass;
           mu1    = body[iBody].dMass / body[0].dMass;
@@ -2557,6 +2557,9 @@ void HessEigen(double **amat, int origsize, double real[], double imag[]) {
                 q /= lrcorner;
                 r /= lrcorner;
               }
+            } else {
+              fprintf(stderr,"ERROR: k=m in distorb.c:HessEigen.");
+              exit(EXIT_INT);
             }
             value = sqrt(p * p + q * q + r * r);
             s     = (double)fiSign(p) * value;
@@ -2653,6 +2656,7 @@ void HessReduce(double **a, int size) {
 
   for (r = 0; r < size; r++) {
     max = 0;
+    rmax = r+1;
     for (rp = r + 1; rp < size; rp++) {
       if (fabs(a[rp][r]) > max) {
         max  = fabs(a[rp][r]);
@@ -3043,6 +3047,9 @@ void RecalcLaplace(BODY *body, EVOLVE *evolve, SYSTEM *system, int iVerbose) {
         alpha1 = body[iBody].dSemi / body[jBody].dSemi;
       } else if (body[iBody].dSemi > body[jBody].dSemi) {
         alpha1 = body[jBody].dSemi / body[iBody].dSemi;
+      } else {
+        fprintf(stderr,"ERROR: Semi-major axes cannot be identical in RecalcLaplace.");
+        exit(EXIT_INPUT);
       }
 
       for (j = 0; j < 26; j++) {
@@ -3076,7 +3083,7 @@ Recalculates eigenvalues in case where LL2 solution is coupled to eqtide
 */
 void RecalcEigenVals(BODY *body, EVOLVE *evolve, SYSTEM *system) {
   int iBody, jBody, j, done = 0;
-  double alpha1, dalpha = -1, dalphaTmp;
+  double alpha1 = 0, dalpha = -1, dalphaTmp;
 
   for (iBody = 1; iBody < evolve->iNumBodies - 1; iBody++) {
     for (jBody = iBody + 1; jBody < evolve->iNumBodies; jBody++) {
@@ -3084,6 +3091,9 @@ void RecalcEigenVals(BODY *body, EVOLVE *evolve, SYSTEM *system) {
         alpha1 = body[iBody].dSemi / body[jBody].dSemi;
       } else if (body[iBody].dSemi > body[jBody].dSemi) {
         alpha1 = body[jBody].dSemi / body[iBody].dSemi;
+      } else {
+        fprintf(stderr,"ERROR: Semi-major axes cannot be identical in RecalcEigenVals.");
+        exit(EXIT_INPUT);
       }
       for (j = 0; j < 2; j++) {
         dalphaTmp =
diff --git a/src/distrot.c b/src/distrot.c
index d9b9c04df..c263d979c 100644
--- a/src/distrot.c
+++ b/src/distrot.c
@@ -419,7 +419,10 @@ void VerifyOrbitData(BODY *body, CONTROL *control, OPTIONS *options,
         exit(EXIT_INPUT);
       }
       // Check file has exactly 7 columns
-      fgets(cLine, LINE, fileorb);
+      if (fgets(cLine, LINE, fileorb) == NULL) {
+        fprintf(stderr,"ERROR: Unable to read line from orbit data file.");
+        exit(EXIT_INPUT);
+      } 
       GetWords(cLine, cFoo, &iNumColsFound, &bFoo);
       if (iNumCols != iNumColsFound) {
         if (control->Io.iVerbose >= VERBERR) {
@@ -456,8 +459,11 @@ void VerifyOrbitData(BODY *body, CONTROL *control, OPTIONS *options,
 
       iLine = 0;
       while (feof(fileorb) == 0) {
-        fscanf(fileorb, "%lf %lf %lf %lf %lf %lf %lf\n", &dttmp, &datmp, &detmp,
-               &ditmp, &daptmp, &dlatmp, &dmatmp);
+        if (fscanf(fileorb, "%lf %lf %lf %lf %lf %lf %lf\n", &dttmp, &datmp, &detmp,
+               &ditmp, &daptmp, &dlatmp, &dmatmp) != 7) {
+                  fprintf(stderr,"ERROR: Incorrect number of columns in orbit file.");
+                  exit(EXIT_INPUT);
+        }
         body[iBody].daTimeSeries[iLine] =
               dttmp * fdUnitsTime(control->Units[iBody + 1].iTime);
         body[iBody].daSemiSeries[iLine] =
diff --git a/src/eqtide.c b/src/eqtide.c
index 04daa4335..9f9bf8e46 100644
--- a/src/eqtide.c
+++ b/src/eqtide.c
@@ -4618,11 +4618,13 @@ double fdDB15DsemiDt(BODY *body, SYSTEM *system, int *iaBody) {
   if (iBody > 0) { // Also means iBody is the orbiter
     double imk2;
     /* Set imk2 using appropriate model */
+    /* XXX This switch is not understood -- deprecate for v3.0 */
     if (body[iBody].dImK2ManOrbModel == 1) {
       imk2 = body[iBody].dImK2Man;
-    }
-    if (body[iBody].dImK2ManOrbModel == 2) {
+    } else if (body[iBody].dImK2ManOrbModel == 2) {
       imk2 = -body[iBody].dK2Man / body[iBody].dTidalQMan;
+    } else {
+      imk2 = body[iBody].dImK2Man; // To eliminate compiler warnings.
     }
     return 21 * imk2 * body[iBody].dSemi * body[iPert].dMass /
            body[iBody].dMass *
@@ -4642,9 +4644,10 @@ double fdDB15DeccDt(BODY *body, UPDATE *update, int *iaBody) {
     /* Set imk2 using appropriate model */
     if (body[iBody].dImK2ManOrbModel == 1) {
       imk2 = body[iBody].dImK2Man;
-    }
-    if (body[iBody].dImK2ManOrbModel == 2) {
+    } else if (body[iBody].dImK2ManOrbModel == 2) {
       imk2 = -body[iBody].dK2Man / body[iBody].dTidalQMan;
+    } else {
+      imk2 = body[iBody].dImK2Man; // To eliminate compiler warnings.
     }
     return 21 / 2. * imk2 * pow(body[iPert].dMass, 3. / 2) * pow(BIGG, 1. / 2) *
            pow(body[iBody].dRadius, 5) * body[iBody].dEcc / body[iBody].dMass /
diff --git a/src/evolve.c b/src/evolve.c
index a93283c14..6aada70bd 100644
--- a/src/evolve.c
+++ b/src/evolve.c
@@ -335,13 +335,9 @@ void fdGetUpdateInfo(BODY *body, CONTROL *control, SYSTEM *system,
 
   int iBody, iVar, iEqn, iNumBodies, iNumVars,
         iNumEqns; // Dummy counting variables
-  EVOLVE
-  integr; // Dummy EVOLVE struct so we don't have to dereference control a lot
   double dVarNow, dMinNow, dMin = dHUGE,
                            dVarTotal; // Intermediate storage variables
 
-  integr = control->Evolve;
-
   iNumBodies = control->Evolve.iNumBodies;
   for (iBody = 0; iBody < iNumBodies; iBody++) {
     if (update[iBody].iNumVars > 0) {
@@ -416,8 +412,6 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
 
   evolve->dCurrentDt = *dDt;
   iNumBodies         = evolve->iNumBodies;
-#pragma omp parallel for num_threads(NUM_THREADS) private(iNumVars, iNumEqns,  \
-                                                          iVar, iEqn)
   for (iBody = 0; iBody < iNumBodies; iBody++) {
     // int thread_num = omp_get_thread_num();
     // int cpu_num = sched_getcpu();
@@ -463,8 +457,6 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
   fdGetUpdateInfo(evolve->tmpBody, control, system, evolve->tmpUpdate,
                   fnUpdate);
 
-#pragma omp parallel for num_threads(NUM_THREADS) private(iNumVars, iNumEqns,  \
-                                                          iVar, iEqn)
   for (iBody = 0; iBody < iNumBodies; iBody++) {
     iNumVars = update[iBody].iNumVars;
     double daDerivVar;
@@ -508,8 +500,6 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
   fdGetUpdateInfo(evolve->tmpBody, control, system, evolve->tmpUpdate,
                   fnUpdate);
 
-#pragma omp parallel for num_threads(NUM_THREADS) private(iNumVars, iNumEqns,  \
-                                                          iVar, iEqn)
   for (iBody = 0; iBody < iNumBodies; iBody++) {
     iNumVars = update[iBody].iNumVars;
     double daDerivVar;
@@ -553,8 +543,6 @@ void RungeKutta4Step(BODY *body, CONTROL *control, SYSTEM *system,
   fdGetUpdateInfo(evolve->tmpBody, control, system, evolve->tmpUpdate,
                   fnUpdate);
 
-#pragma omp parallel for num_threads(NUM_THREADS) private(iNumVars, iNumEqns,  \
-                                                          iVar, iEqn)
   for (iBody = 0; iBody < iNumBodies; iBody++) {
     double daDerivVar;
     iNumVars = update[iBody].iNumVars;
diff --git a/src/flare.c b/src/flare.c
index 1e28291e0..3489fd84b 100644
--- a/src/flare.c
+++ b/src/flare.c
@@ -755,7 +755,7 @@ void VerifyFlare(BODY *body,
   already inputed the XUV luminosity by flares and the flare module doesn't need
   this information anymore*/
   if (body[iBody].iFlareFFD == FLARE_FFD_NONE) {
-    int iCol, bError = 0;
+    int iCol;
     for (iCol = 0; iCol < files->Outfile[iBody].iNumCols; iCol++) {
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREFREQ1].cName,
@@ -765,7 +765,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREFREQ1].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREFREQ2].cName,
@@ -775,7 +774,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREFREQ2].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREFREQ3].cName,
@@ -785,7 +783,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREFREQ3].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREFREQ4].cName,
@@ -795,7 +792,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREFREQ4].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREFREQMIN].cName,
@@ -805,7 +801,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREFREQMIN].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREFREQMID].cName,
@@ -815,7 +810,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREFREQMID].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREFREQMAX].cName,
@@ -825,7 +819,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREFREQMAX].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREENERGY1].cName,
@@ -835,7 +828,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREENERGY1].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREENERGY2].cName,
@@ -845,7 +837,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREENERGY2].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREENERGY3].cName,
@@ -855,7 +846,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREENERGY3].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREENERGY4].cName,
@@ -865,7 +855,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREENERGY4].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREENERGYMIN].cName,
@@ -875,7 +864,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREENERGYMIN].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREENERGYMID].cName,
@@ -885,7 +873,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREENERGYMID].cName);
-        bError = 1;
       }
       if (memcmp(files->Outfile[iBody].caCol[iCol],
                  output[OUT_FLAREENERGYMAX].cName,
@@ -895,7 +882,6 @@ void VerifyFlare(BODY *body,
                 "WARNING: Output option %s only allowed with FFD model "
                 "DAVENPORT or LACY \n",
                 output[OUT_FLAREENERGYMAX].cName);
-        bError = 1;
       }
     }
   }
diff --git a/src/galhabit.c b/src/galhabit.c
index e5908b764..e08a7b26b 100644
--- a/src/galhabit.c
+++ b/src/galhabit.c
@@ -1845,18 +1845,13 @@ void AddModuleGalHabit(CONTROL *control, MODULE *module, int iBody,
 /************* GALHABIT Functions ***********/
 void PropsAuxGalHabit(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update,
                       int iBody) {
-  double sinw, cosw, cosw_alt, sign, dMu, dL;
+  double sinw, cosw, cosw_alt, sign;
 
   /* calculate osculating elements */
   body[iBody].dEcc =
         sqrt(pow(body[iBody].dEccX, 2) + pow(body[iBody].dEccY, 2) +
              pow(body[iBody].dEccZ, 2));
-  //   body[iBody].dAngM = sqrt(pow(body[iBody].dAngMX,2)+pow(body[iBody].dAngMY,2)+\
-                                                          pow(body[iBody].dAngMZ,2));
-  dMu = BIGG * (body[iBody].dMassInterior +
-                body[iBody].dMass); // calculate mass coefficient for
-                                    // primary/primary+secondary
-  dL                = sqrt(dMu * body[iBody].dSemi);
+ 
   body[iBody].dAngM = sqrt(1.0 - pow(body[iBody].dEcc, 2));
 
   body[iBody].dInc   = acos(body[iBody].dAngMZ / body[iBody].dAngM);
@@ -1875,10 +1870,6 @@ void PropsAuxGalHabit(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update,
 
   body[iBody].dArgP = atan2(sinw, cosw);
 
-  // if (body[iBody].bHostBinary) {
-  //     Rot2Bin(body,iBody);
-  //   }
-
   body[iBody].dLongP = body[iBody].dLongA + body[iBody].dArgP;
   while (body[iBody].dArgP > 2 * PI) {
     body[iBody].dArgP -= 2 * PI;
@@ -1898,8 +1889,6 @@ void PropsAuxGalHabit(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update,
   while (body[iBody].dLongA < 0) {
     body[iBody].dLongA += 2 * PI;
   }
-
-  //   body[iBody].dEcc = 1.0 - body[iBody].dPeriQ/body[iBody].dSemi;
 }
 
 void ForceBehaviorGalHabit(BODY *body, MODULE *module, EVOLVE *evolve, IO *io,
@@ -2256,10 +2245,6 @@ void CalcEccVec(BODY *body, int iBody) {
 }
 
 void CalcAngMVec(BODY *body, int iBody) {
-  double dMu, dL;
-  dMu = BIGG * (body[iBody].dMassInterior); //+body[iBody].dMass);
-  dL  = sqrt(dMu * body[iBody].dSemi);
-
   body[iBody].dAngM = sqrt((1.0 - pow(body[iBody].dEcc, 2)));
   body[iBody].dAngMX =
         body[iBody].dAngM * (sin(body[iBody].dLongA) * sin(body[iBody].dInc));
@@ -2998,9 +2983,7 @@ double fndGalHabitDEccZDtTidal(BODY *body, SYSTEM *system, int *iaBody) {
 }
 
 double fndGalHabitDAngMXDtTidal(BODY *body, SYSTEM *system, int *iaBody) {
-  double dMu, dJ, dL;
-  dMu = BIGG * (body[iaBody[0]].dMassInterior); //+body[iaBody[0]].dMass);
-  dL  = sqrt(dMu * body[iaBody[0]].dSemi);
+  double dJ;
   dJ  = sqrt((1.0 - pow(body[iaBody[0]].dEcc, 2)));
 
   return sin(body[iaBody[0]].dLongA) * sin(body[iaBody[0]].dInc) *
@@ -3010,9 +2993,7 @@ double fndGalHabitDAngMXDtTidal(BODY *body, SYSTEM *system, int *iaBody) {
 }
 
 double fndGalHabitDAngMYDtTidal(BODY *body, SYSTEM *system, int *iaBody) {
-  double dMu, dJ, dL;
-  dMu = BIGG * (body[iaBody[0]].dMassInterior); //+body[iaBody[0]].dMass);
-  dL  = sqrt(dMu * body[iaBody[0]].dSemi);
+  double dJ;
   dJ  = sqrt((1.0 - pow(body[iaBody[0]].dEcc, 2)));
 
   return -cos(body[iaBody[0]].dLongA) * sin(body[iaBody[0]].dInc) *
diff --git a/src/magmoc.c b/src/magmoc.c
index e54108125..562a8a3b1 100644
--- a/src/magmoc.c
+++ b/src/magmoc.c
@@ -845,25 +845,30 @@ double fndBisection(double (*f)(BODY *, double, int), BODY *body, double dXl,
                     double dXu, double dEps, int iBody) {
   double dXm, dEpsilon, dProd, dFxm, dFxl;
   dEpsilon = 10 * dEps;
-  while (dEpsilon > dEps) {
-    dXm  = (dXl + dXu) / 2.;
-    dFxm = (*f)(body, dXm, iBody);
-    if (fabs(dFxm) < dEps) {
-      return dXm;
-    }
-    dFxl = (*f)(body, dXl, iBody);
-    if (fabs(dFxl) < dEps) {
-      return dXl;
-    }
-    dProd = (dFxl / fabs(dFxl)) * (dFxm / fabs(dFxm));
-    if (dProd < 0) {
-      dXu = dXm;
-    } else {
-      dXl = dXm;
+  if (dEpsilon > dEps) {
+    while (dEpsilon > dEps) {
+      dXm  = (dXl + dXu) / 2.;
+      dFxm = (*f)(body, dXm, iBody);
+      if (fabs(dFxm) < dEps) {
+        return dXm;
+      }
+      dFxl = (*f)(body, dXl, iBody);
+      if (fabs(dFxl) < dEps) {
+        return dXl;
+      }
+      dProd = (dFxl / fabs(dFxl)) * (dFxm / fabs(dFxm));
+      if (dProd < 0) {
+        dXu = dXm;
+      } else {
+        dXl = dXm;
+      }
+      dEpsilon = fabs((*f)(body, dXm, iBody));
     }
-    dEpsilon = fabs((*f)(body, dXm, iBody));
-  }
-  return dXm;
+    return dXm;
+  } else {
+    fprintf(stderr,"ERROR: Tolerance factor <= 0 in fndBisection.");
+    exit(EXIT_INT);
+  }  
 }
 
 /**
diff --git a/src/options.c b/src/options.c
index da151ceed..2e6557aae 100644
--- a/src/options.c
+++ b/src/options.c
@@ -100,7 +100,10 @@ void GetNextValidLine(char cFile[], int iStart, char cLine[], int *iLine) {
   *iLine = 0;
 
   for (iLineTmp = 0; iLineTmp < iStart; iLineTmp++) {
-    fgets(cLine, LINE, fp);
+    if (fgets(cLine, LINE, fp) == NULL) {
+      fprintf(stderr,"ERROR: Unable to read next valid line.");
+      LineExit(cFile,iStart); 
+    }
     (*iLine)++;
   }
 
@@ -373,7 +376,7 @@ void AddOptionString(char cFile[], char cOption[], char cInput[], int *iLine,
   sscanf(cLine, "%s %s", cTmp, cInput);
 }
 
-
+/* Looks like this was deprecated somewhere RB 01/02/24 
 int GetNumOut(char cFile[], char cName[], int iLen, int *iLineNum, int iExit) {
   char cLine[LINE], cWord[NAMELEN];
   int iPos, j, ok, bDone = 0, iLine = 0, iNumOut;
@@ -386,12 +389,12 @@ int GetNumOut(char cFile[], char cName[], int iLen, int *iLineNum, int iExit) {
   }
 
   while (fgets(cLine, LINE, fp) != NULL) {
-    /* Check for # sign */
+    // Check for # sign 
     if (memcmp(cLine, "#", 1) != 0) {
-      /* Check for desired parameter */
+      // Check for desired parameter 
       sscanf(cLine, "%s", cWord);
       if (memcmp(cWord, cName, iLen) == 0) {
-        /* Parameter Found! */
+        // Parameter Found! 
         if (bDone) {
           fprintf(stderr, "ERROR: Multiple occurences of parameter %s found.\n",
                   cName);
@@ -404,12 +407,12 @@ int GetNumOut(char cFile[], char cName[], int iLen, int *iLineNum, int iExit) {
         iNumOut = 0;
         ok      = 1;
         for (iPos = 1; iPos < LINE;
-             iPos++) { /* Ignore first character, as it makes conditional
-                          well-defined */
-          /* printf("%d ",cLine[iPos]); */
+             iPos++) { // Ignore first character, as it makes conditional
+                          // well-defined 
+          // printf("%d ",cLine[iPos]); 
           if (ok) {
             if (cLine[iPos] == 35) { // 35 is ASCII code for #
-              /* Pound sign! */
+              // Pound sign! 
               ok = 0;
               iNumOut++;
             }
@@ -425,10 +428,11 @@ int GetNumOut(char cFile[], char cName[], int iLen, int *iLineNum, int iExit) {
       cLine[iPos] = 0;
     }
   }
-  /* Lose the input parameter */
+  // Lose the input parameter
   iNumOut--;
   return iNumOut;
 }
+*/
 
 int iGetNumLines(char cFile[]) {
   int iNumLines = 0, iChar, bFileOK = 1;
@@ -487,7 +491,7 @@ void InitializeInput(INFILE *input) {
 
   fp = fopen(input->cIn, "r");
   if (fp == NULL) {
-    fprintf(stderr, "Unable to open %s.\n", input->cIn);
+    fprintf(stderr, "ERROR: Unable to open %s.\n", input->cIn);
     exit(EXIT_INPUT);
   }
   input->iNumLines = iGetNumLines(input->cIn);
@@ -501,7 +505,10 @@ void InitializeInput(INFILE *input) {
     input->bLineOK[iLine] = 0;
     memset(cLine, '\0', LINE);
 
-    fgets(cLine, LINE, fp);
+    if (fgets(cLine, LINE, fp) == NULL) {
+      fprintf(stderr,"ERROR: Unable to open %s.\n", input->cIn);
+      exit(EXIT_INPUT);
+    }
     if (CheckComment(cLine, LINE)) {
       input->bLineOK[iLine] = 1;
     } else {
@@ -2643,7 +2650,7 @@ void ReadCosObl(BODY *body, CONTROL *control, FILES *files, OPTIONS *options,
 void ReadOutputOrder(FILES *files, MODULE *module, OPTIONS *options,
                      OUTPUT *output, int iFile, int iVerbose) {
   int i, j, count, iLen, iNumIndices = 0, bNeg[MAXARRAY], ok = 1, iNumGrid = 0;
-  int k, iOut, *lTmp, iCol, jCol;
+  int k, iOut=-1, *lTmp, iCol, jCol;
   char saTmp[MAXARRAY][OPTLEN], cTmp[OPTLEN], cOption[MAXARRAY][OPTLEN],
         cOut[OPTLEN];
   int iLen1, iLen2;
@@ -2846,7 +2853,7 @@ void ReadOutputOrder(FILES *files, MODULE *module, OPTIONS *options,
 void ReadGridOutput(FILES *files, OPTIONS *options, OUTPUT *output, int iFile,
                     int iVerbose) {
   int i, j, count, iLen, iNumIndices = 0, bNeg[MAXARRAY], ok = 0, iNumGrid = 0;
-  int k, iOut, *lTmp;
+  int k, iOut=-1, *lTmp;
   char saTmp[MAXARRAY][OPTLEN], cTmp[OPTLEN], cOption[MAXARRAY][OPTLEN],
         cOut[OPTLEN];
   int iLen1, iLen2;
diff --git a/src/output.c b/src/output.c
index 487df799e..3b6bb76cd 100644
--- a/src/output.c
+++ b/src/output.c
@@ -2003,9 +2003,9 @@ void LogOutputOrder(BODY *body, CONTROL *control, FILES *files, OUTPUT *output,
                     SYSTEM *system, UPDATE *update, fnWriteOutput fnWrite[],
                     FILE *fp, int iBody) {
   int iCol, iOut, iSubOut, iExtra = 0;
-  char cCol[MODULEOUTEND][OUTLEN];
+  char cCol[MODULEOUTEND][OUTLEN+2]; // +2 for brackets
   double *dTmp;
-  char cUnit[OUTLEN], cTmp[OUTLEN];
+  char cUnit[OUTLEN], cTmp[OUTLEN+2];
 
   for (iCol = 0; iCol < files->Outfile[iBody].iNumCols; iCol++) {
     for (iOut = 0; iOut < MODULEOUTEND; iOut++) {
@@ -2038,9 +2038,9 @@ void LogGridOutput(BODY *body, CONTROL *control, FILES *files, OUTPUT *output,
                    SYSTEM *system, UPDATE *update, fnWriteOutput fnWrite[],
                    FILE *fp, int iBody) {
   int iCol, iOut, iSubOut, iExtra = 0;
-  char cCol[MODULEOUTEND][OUTLEN];
+  char cCol[MODULEOUTEND][OUTLEN+2]; // +2 for brackets
   double *dTmp;
-  char cUnit[OUTLEN], cTmp[OUTLEN];
+  char cUnit[OUTLEN], cTmp[OUTLEN+2];
 
   for (iCol = 0; iCol < files->Outfile[iBody].iNumGrid; iCol++) {
     for (iOut = 0; iOut < MODULEOUTEND; iOut++) {
@@ -2163,11 +2163,14 @@ void WriteLog(BODY *body, CONTROL *control, FILES *files, MODULE *module,
               fnUpdateVariable ***fnUpdate, fnWriteOutput fnWrite[], int iEnd) {
   char cTime[OPTLEN];
   FILE *fp;
-  double dDt, dTotTime;
+  double dTotTime;
 
   /* Get derivatives */
   PropertiesAuxiliary(body, control, system, update);
-  dDt = fdGetTimeStep(body, control, system, update, fnUpdate);
+  /* XXX The fdGetTimeStep function is not single-purpose. The function sets 
+     members of the UPDATE struct, so those value should probably  be set in 
+     one function, and then the timestep calculated in another. */
+  double dDt = fdGetTimeStep(body, control, system, update, fnUpdate);
 
   if (iEnd == 0) {
     sprintf(cTime, "Input");
diff --git a/src/poise.c b/src/poise.c
index f0cc4e5da..d3cfc28f3 100644
--- a/src/poise.c
+++ b/src/poise.c
@@ -1824,8 +1824,11 @@ void VerifyOrbitOblData(BODY *body, CONTROL *control, OPTIONS *options,
 
       iLine = 0;
       while (feof(fileorb) == 0) {
-        fscanf(fileorb, "%lf %lf %lf %lf %lf %lf %lf", &dttmp, &datmp, &detmp,
-               &daptmp, &dlatmp, &dobltmp, &dprecatmp);
+        if (fscanf(fileorb, "%lf %lf %lf %lf %lf %lf %lf", &dttmp, &datmp, &detmp,
+              &daptmp, &dlatmp, &dobltmp, &dprecatmp) != 7) {
+                fprintf(stderr,"ERROR: Incorrect number of columns in orbit-obliquity file.");
+                exit(EXIT_INPUT);
+        }
 
         body[iBody].daTimeSeries[iLine] =
               dttmp * fdUnitsTime(control->Units[iBody + 1].iTime);
@@ -1964,7 +1967,7 @@ void DampTemp(BODY *body, double dTGlobalTmp, int iBody) {
 void InitializeClimateParams(BODY *body, int iBody, int iVerbose) {
   int iLat, jLat, count, iRun;
   double Toffset, dXBoundary, RunningMeanTmp;
-  double *daRunningMean, TotalMean;
+  double *daRunningMean;
 
   body[iBody].dIceMassTot   = 0.0;
   body[iBody].daInsol       = malloc(body[iBody].iNumLats * sizeof(double *));
@@ -2451,7 +2454,6 @@ void InitializeClimateParams(BODY *body, int iBody, int iVerbose) {
       int iMaxIteration     = 2 * RunLen;
       daRunningMean         = malloc((RunLen + 1) * sizeof(double));
       daRunningMean[RunLen] = 0;
-      TotalMean             = 0;
       RunningMeanTmp        = daRunningMean[RunLen];
       while (fabs(RunningMeanTmp - daRunningMean[RunLen]) >
                    body[iBody].dSpinUpTol ||
@@ -2623,10 +2625,6 @@ void NullPoiseDerivatives(BODY *body, EVOLVE *evolve, UPDATE *update,
 void VerifyPoise(BODY *body, CONTROL *control, FILES *files, OPTIONS *options,
                  OUTPUT *output, SYSTEM *system, UPDATE *update, int iBody,
                  int iModule) {
-  int iLat, jLat;
-  jLat = 0;
-  iLat = 0;
-
   VerifyAlbedo(body, options, files->Infile[iBody + 1].cIn, iBody,
                control->Io.iVerbose);
   VerifyAstro(body, options, files->Infile[iBody + 1].cIn, iBody,
@@ -4789,8 +4787,7 @@ cap.
 void fvNorthIceCapLand(BODY *body, int iBody, double *dLatIceEdge,
                        int *iLatIceEdge, int *bCap) {
 
-  int iLat, iNum;
-  iNum = 0;
+  int iLat;
 
   // Check for ice at north pole; no ice at +90 => No ice cap
   if (!fbIceLatLand(body, iBody, body[iBody].iNumLats - 1)) {
@@ -4835,8 +4832,7 @@ Determines if planet has a northern polar sea ice cap and the extent of the cap.
 void fvNorthIceCapSea(BODY *body, int iBody, double *dLatIceEdge,
                       int *iLatIceEdge, int *bCap) {
 
-  int iLat, iNum;
-  iNum = 0;
+  int iLat;
 
   // Check for ice at north pole; no ice at +90 => No ice cap
   if (!fbIceLatSea(body, iBody, body[iBody].iNumLats - 1)) {
@@ -4882,8 +4878,7 @@ cap.
 void fvSouthIceCapLand(BODY *body, int iBody, double *dLatIceEdge,
                        int *iLatIceEdge, int *bCap) {
 
-  int iLat, iNum;
-  iNum = 0;
+  int iLat;
 
   // Check for ice at south pole; no ice at -90 => No ice cap
   if (!fbIceLatLand(body, iBody, 0)) {
@@ -4928,8 +4923,7 @@ Determines if planet has a southern polar sea ice cap
 void fvSouthIceCapSea(BODY *body, int iBody, double *dLatIceEdge,
                       int *iLatIceEdge, int *bCap) {
 
-  int iLat, iNum;
-  iNum = 0;
+  int iLat;
 
   // Check for ice at south pole; no ice at -90 => No ice cap
   if (!fbIceLatSea(body, iBody, 0)) {
@@ -4978,7 +4972,7 @@ void fvIceBeltLand(BODY *body, int iBody, double *dLatIceEdgeNorth,
                    int *iLatIceEdgeSouth, int *bBelt) {
 
   int bCapNorth, bCapSouth, iIce;
-  int iLat, iEquator, iLatStart, iLatEnd, iLatCapNorth, iLatCapSouth;
+  int iLat, iLatStart, iLatEnd, iLatCapNorth, iLatCapSouth;
   double dLatCapNorth, dLatCapSouth;
 
   // If IceFree or Snowball, no ice belt
@@ -4998,7 +4992,6 @@ void fvIceBeltLand(BODY *body, int iBody, double *dLatIceEdgeNorth,
   // If made it here, belt is possible, so let's look!
   *iLatIceEdgeNorth = 0;
   *iLatIceEdgeSouth = 0;
-  iEquator          = (int)(body[iBody].iNumLats / 2);
 
   // Get parameters for Northern Ice Cap
   fvNorthIceCapLand(body, iBody, &dLatCapNorth, &iLatCapNorth, &bCapNorth);
@@ -5738,10 +5731,9 @@ operation of the inverse matrix is a "time-step".
 */
 void fvMatrixAnnual(BODY *body, int iBody) {
   int iLat, jLat;
-  double dDelta_t, dDelta_x;
+  double dDelta_t;
 
   dDelta_t = 1.5 / body[iBody].iNumLats;
-  dDelta_x = 2.0 / body[iBody].iNumLats;
 
   for (iLat = 0; iLat < body[iBody].iNumLats; iLat++) {
     body[iBody].daTempTerms[iLat] = 0.0;
@@ -6748,14 +6740,13 @@ operation of the inverse matrix is a "time-step".
 */
 void fvMatrixSeasonal(BODY *body, int iBody) {
   int iLat, jLat;
-  double dXBoundary, dNu_fl, dNu_fw, dCl_dt, dCw_dt, dt_Lf;
+  double dXBoundary, dNu_fl, dNu_fw, dCl_dt, dCw_dt;
 
   dCl_dt = body[iBody].dHeatCapLand * body[iBody].dMeanMotion / (2 * PI) /
            body[iBody].dSeasDeltat;
 
   dCw_dt = body[iBody].dHeatCapWater * body[iBody].dMeanMotion / (2 * PI) /
            body[iBody].dSeasDeltat;
-  dt_Lf = body[iBody].dSeasDeltat / body[iBody].dLatentHeatIce;
 
   for (iLat = 0; iLat < body[iBody].iNumLats + 1; iLat++) {
     body[iBody].daLambdaSea[iLat] =
@@ -6873,15 +6864,13 @@ should be used as a check on the physics in the EBM (should be ~0).
 */
 void EnergyResiduals(BODY *body, int iBody, int iNday) {
   int iLat;
-  double dNu_fl, dNu_fw, dCl_dt, dCw_dt, dt_Lf;
+  double dNu_fl, dNu_fw, dCl_dt, dCw_dt;
 
   dCl_dt = body[iBody].dHeatCapLand * body[iBody].dMeanMotion / (2 * PI) /
            body[iBody].dSeasDeltat;
 
   dCw_dt = body[iBody].dHeatCapWater * body[iBody].dMeanMotion / (2 * PI) /
            body[iBody].dSeasDeltat;
-  dt_Lf = body[iBody].dSeasDeltat * body[iBody].dMeanMotion / (2 * PI) /
-          body[iBody].dLatentHeatIce;
 
   for (iLat = 0; iLat < body[iBody].iNumLats; iLat++) {
     dNu_fl = body[iBody].dNuLandWater / body[iBody].daLandFrac[iLat];
diff --git a/src/radheat.c b/src/radheat.c
index fae554650..0145786dd 100644
--- a/src/radheat.c
+++ b/src/radheat.c
@@ -3588,8 +3588,6 @@ int fbHaltMin235UPower(BODY *body, EVOLVE *evolve, HALT *halt, IO *io,
 */
 int fbHaltMinRadPower(BODY *body, EVOLVE *evolve, HALT *halt, IO *io,
                       UPDATE *update, fnUpdateVariable ***fnUpdate, int iBody) {
-  int iFoo;
-  iFoo = fdRadPowerTotal(body, iBody);
 
   if (fdRadPowerTotal(body, iBody) < halt->dMinRadPower) {
     if (io->iVerbose >= VERBPROG) {
diff --git a/src/stellar.c b/src/stellar.c
index 91df37b1a..f00d3d40b 100644
--- a/src/stellar.c
+++ b/src/stellar.c
@@ -757,7 +757,7 @@ void fnPropsAuxStellar(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update,
   if (body[iBody].iXUVModel == STELLAR_MODEL_REINERS) {
 
     // REINERS wind model
-    double dPer, dLXRay, dLXRaySat, dLEUV;
+    double dPer, dLXRay, dLXRaySat;
     dPer = 2 * PI / body[iBody].dRotRate;
 
     // Unsaturated regime (Reiners, Schussler & Passegger 2014, eqn. (11))
@@ -772,8 +772,10 @@ void fnPropsAuxStellar(BODY *body, EVOLVE *evolve, IO *io, UPDATE *update,
       dLXRay = dLXRaySat;
     }
 
-    // Sanz-Forcada et al. (2011), eqn (3)
+    /* Sanz-Forcada et al. (2011), eqn (3)
+    Not used here, but maybe useful elsewhere?
     dLEUV = 1.e7 * pow(10., 4.80 + 0.860 * log10(dLXRay * 1.e-7));
+    */
 
     // NOTE: We should add XRay and EUV to get XUV, but the Sanz-Forcada
     // model above yields unrealistically high EUV luminosities for M dwarfs.
diff --git a/src/vplanet.h b/src/vplanet.h
index 012d6d689..06037e0a5 100644
--- a/src/vplanet.h
+++ b/src/vplanet.h
@@ -1911,8 +1911,8 @@ struct INFILE {
  * regarding the output files. */
 
 struct OUTFILE {
-  char cOut[NAMELEN];                /**< Output File Name */
-  int iNumCols;                      /**< Number of Columns in Output File */
+  char cOut[2*NAMELEN+10];           /**< Output File Name */
+  int iNumCols;                      /**< Number of Columns in Output File (system.planet+.forward/backward) */
   char caCol[MODULEOUTEND][OPTLEN];  /**< Output Value Name */
   int bNeg[MODULEOUTEND];            /**< Use Negative Option Units? */
   int iNumGrid;                      /**< Number of grid outputs */
@@ -1926,7 +1926,7 @@ struct OUTFILE {
 struct FILES {
   char cExe[LINE];    /**< Name of Executable */
   OUTFILE *Outfile;   /**< Output File Name for Forward Integration */
-  char cLog[NAMELEN]; /**< Log File Name */
+  char cLog[NAMELEN+4]; /**< Log File Name (+4 to allow for ".log" suffix) */
   INFILE *Infile;
   int iNumInputs; /**< Number of Input Files */
 };