Implemented the KCL verification in place of the classical SPICE3 check for the circuit currents

At the moment the used formula is not the best approach, but it avoids to rewrite every model

src/include/ngspice/cktdefs.h

@@ -425,6 +425,10 @@ extern void NInzIter(CKTcircuit *, int, int);
 extern int NIpred(CKTcircuit *ckt);
 #endif
 
+/* KCL Verification */
+extern int NIkclVerification (CKTcircuit *) ;
+/********************/
+
 extern IFfrontEnd *SPfrontEnd;
 
 #endif

src/include/ngspice/smpdefs.h

@@ -14,6 +14,10 @@ Modified: 2000  AlansFixes
 #include <math.h>
 #include "ngspice/complex.h"
 
+/* KCL Verification */
+int KCLverification (SMPmatrix *, double *, double *, double, int, double, double, double) ;
+/********************/
+
 int SMPaddElt( SMPmatrix *, int , int , double );
 double * SMPmakeElt( SMPmatrix * , int , int );
 void SMPcClear( SMPmatrix *);

src/include/ngspice/spmatrix.h

@@ -294,4 +294,8 @@ extern  void     spMultTransposed(MatrixPtr,spREAL*,spREAL*,spREAL*,spREAL*);
 extern  void     spSolve( MatrixPtr, spREAL*, spREAL*, spREAL*, spREAL* );
 extern  void     spSolveTransposed(MatrixPtr,spREAL*,spREAL*,spREAL*,spREAL*);
 
+/* KCL Verification */
+extern int KCL_verification (MatrixPtr, spREAL *, spREAL *, double, int, double, double, double) ;
+/********************/
+
 #endif  /* spOKAY */

src/maths/ni/niconv.c

@@ -45,18 +45,6 @@ NIconvTest(CKTcircuit *ckt)
 #ifdef STEPDEBUG
                 printf(" non-convergence at node (type=3) %s (fabs(new-old)>tol --> fabs(%g-%g)>%g)\n",CKTnodName(ckt,i),new,old,tol);
 		printf("    reltol: %g    voltTol: %g   (tol=reltol*(MAX(fabs(old),fabs(new))) + voltTol)\n",ckt->CKTreltol,ckt->CKTvoltTol);
-#endif /* STEPDEBUG */
-		ckt->CKTtroubleNode = i;
-		ckt->CKTtroubleElt = NULL;
-                return(1);
-            }
-        } else {
-            tol =  ckt->CKTreltol * (MAX(fabs(old),fabs(new))) +
-                    ckt->CKTabstol;
-            if (fabs(new-old) >tol ) {
-#ifdef STEPDEBUG
-                printf(" non-convergence at node (type=%d) %s (fabs(new-old)>tol --> fabs(%g-%g)>%g)\n",node->type,CKTnodName(ckt,i),new,old,tol);
-		printf("    reltol: %g    abstol: %g   (tol=reltol*(MAX(fabs(old),fabs(new))) + abstol)\n",ckt->CKTreltol,ckt->CKTabstol);
 #endif /* STEPDEBUG */
 		ckt->CKTtroubleNode = i;
 		ckt->CKTtroubleElt = NULL;
@@ -74,3 +62,41 @@ NIconvTest(CKTcircuit *ckt)
     return(0);
 #endif /* NEWCONV */
 }
+
+
+/**
+ * Routine to Verify the KCL
+ */
+
+int NIkclVerification (CKTcircuit *ckt)
+{
+    int i, size ;
+    double maximum = 0 ;
+    CKTnode *node ;
+
+    size = SMPmatSize (ckt->CKTmatrix) ;
+
+    node = ckt->CKTnodes ;
+    for (i = 1 ; i <= size ; i++)
+    {
+        if (node->type == SP_CURRENT)
+        {
+            if (maximum < fabs (ckt->CKTrhsOld [i]))
+                maximum = fabs (ckt->CKTrhsOld [i]) ;
+        }
+        node = node->next ;
+    }
+
+    node = ckt->CKTnodes ;
+    for (i = 1 ; i <= size ; i++)
+    {
+        if (node->type == SP_VOLTAGE)
+        {
+            if (KCLverification (ckt->CKTmatrix, ckt->CKTrhsOld, ckt->CKTrhs, ckt->CKTdiagGmin, i, ckt->CKTreltol, ckt->CKTabstol, maximum))
+                return 1 ;
+        }
+        node = node->next ;
+    }
+
+    return 0 ;
+}

src/maths/ni/niiter.c

@@ -37,6 +37,9 @@ NIiter(CKTcircuit *ckt, int maxIter)
     CKTnode *node; /* current matrix entry */
     double diff, maxdiff, damp_factor, *OldCKTstate0=NULL;
 
+    /* KCL_verified */
+    int KCL_verified ;
+
     if ( maxIter < 100 ) maxIter = 100; /* some convergence issues that get resolved by increasing max iter */
 
     iterno=0;
@@ -110,6 +113,16 @@ NIiter(CKTcircuit *ckt, int maxIter)
                 ckt->CKTniState |= NISHOULDREORDER;
             }
 
+            /* KCL Verification */
+            if ((ckt->CKTnoncon == 0) && (iterno != 1))
+            {
+                KCL_verified = NIkclVerification (ckt) ;
+                if (KCL_verified) /* NOT VERIFIED */
+                    ckt->CKTnoncon = 1 ;
+                else
+                    ckt->CKTnoncon = 0 ;
+            }
+
             if(ckt->CKTniState & NISHOULDREORDER) {
                 startTime = SPfrontEnd->IFseconds();
                 error = SMPreorder(ckt->CKTmatrix,ckt->CKTpivotAbsTol,

src/maths/sparse/Makefile.am

@@ -9,6 +9,7 @@ libsparse_la_SOURCES = \
 	spdefs.h	\
 	spextra.c	\
 	spfactor.c	\
+	spkirch.c	\
 	spoutput.c	\
 	spsmp.c		\
 	spsolve.c	\

src/maths/sparse/spkirch.c

@@ -0,0 +1,41 @@
+/**
+ * Routine to Verify the KCL
+ */
+
+#include <math.h>
+#define spINSIDE_SPARSE
+#include "spconfig.h"
+#include "ngspice/spmatrix.h"
+#include "spdefs.h"
+
+int KCL_verification (MatrixPtr Matrix, spREAL *rhsOld, spREAL *rhs, double gmin, int i, double RelTol, double AbsTol, double maximum)
+{
+    ElementPtr element ;
+    spREAL current ;
+
+    current = 0 ;
+
+#ifdef TRANSLATE
+    element = Matrix->FirstInRow [Matrix->ExtToIntRowMap [i]] ;
+#else
+    element = Matrix->FirstInRow [i] ;
+#endif
+
+    /* A*x */
+    while (element != NULL)
+    {
+        current += element->Real * rhsOld [Matrix->IntToExtColMap [element->Col]] ;
+        element = element->NextInRow ;
+    }
+
+#ifdef TRANSLATE
+    current += gmin * rhsOld [Matrix->IntToExtColMap [Matrix->Diag [Matrix->ExtToIntRowMap [i]]->Col]] ;
+#else
+    current += gmin * rhsOld [Matrix->IntToExtColMap [Matrix->Diag [i]->Col]] ;
+#endif
+
+    if (fabs (current - rhs [i]) > RelTol * maximum + AbsTol)
+        return 1 ;
+
+    return 0 ;
+}

src/maths/sparse/spsmp.c

@@ -108,6 +108,14 @@ extern double logb(double);
 
 static void LoadGmin(SMPmatrix *eMatrix, double Gmin);
 
+/**
+ * Routine to Verify the KCL
+ */
+int
+KCLverification (SMPmatrix *Matrix, double *rhsOld, double *rhs, double gmin, int i, double RelTol, double AbsTol, double maximum)
+{
+    return KCL_verification (Matrix, rhsOld, rhs, gmin, i, RelTol, AbsTol, maximum) ;
+}
 
 /*
  * SMPaddElt()