diff --git a/src/frontend/commands.c b/src/frontend/commands.c
index 548c5459c..cbc7c1f5a 100644
--- a/src/frontend/commands.c
+++ b/src/frontend/commands.c
@@ -293,6 +293,14 @@ struct comm spcp_coms[] = {
       "[.pss line args] : Do a periodic state analysis." } ,
 /* SP */
 #endif
+
+//#ifdef RELAN
+    { "relan", com_relan, TRUE, TRUE,
+      { 0, 0, 0, 0 }, E_DEFHMASK, 0, LOTS,
+      NULL,
+      "[.relan line args] : Do a Reliability Analysis." } ,
+//#endif
+
     { "ac", com_ac, TRUE, TRUE,
       { 0, 0, 0, 0 }, E_DEFHMASK, 0, LOTS,
       NULL,
@@ -724,6 +732,14 @@ struct comm nutcp_coms[] = {
       "[.pss line args] : Do a periodic steady state analysis." } ,
 /* SP */
 #endif
+
+//#ifdef RELAN
+    { "relan", NULL, TRUE, TRUE,
+      { 0, 0, 0, 0 }, E_DEFHMASK, 0, LOTS,
+      NULL,
+      "[.relan line args] : Do a Reliability Analysis." } ,
+//#endif
+
     { "ac", NULL, TRUE, TRUE,
       { 0, 0, 0, 0 }, E_DEFHMASK, 0, LOTS,
       NULL,
diff --git a/src/frontend/runcoms.c b/src/frontend/runcoms.c
index 2fd724083..f9f8bfb74 100644
--- a/src/frontend/runcoms.c
+++ b/src/frontend/runcoms.c
@@ -185,6 +185,15 @@ com_pss(wordlist *wl)
 #endif
 
 
+//#ifdef RELAN
+void
+com_relan (wordlist *wl)
+{
+    dosim ("relan", wl) ;
+}
+//#endif
+
+
 static int
 dosim(
     char *what, /* in: command (pz,op,dc,ac,tf,tran,sens,disto,noise,run) */
diff --git a/src/frontend/runcoms.h b/src/frontend/runcoms.h
index 714239e50..364a66fd8 100644
--- a/src/frontend/runcoms.h
+++ b/src/frontend/runcoms.h
@@ -16,6 +16,11 @@ void com_tran(wordlist *wl);
 /* SP: Stady State Analysis */
 void com_pss(wordlist *wl);
 /* SP */
+
+//#ifdef RELAN
+void com_relan (wordlist *wl) ;
+//#endif
+
 void com_sens(wordlist *wl);
 void com_disto(wordlist *wl);
 void com_noise(wordlist *wl);
diff --git a/src/frontend/shyu.c b/src/frontend/shyu.c
index ef2f65a9f..1c5bc02f2 100644
--- a/src/frontend/shyu.c
+++ b/src/frontend/shyu.c
@@ -314,6 +314,70 @@ if_sens_run(CKTcircuit *ckt, wordlist *args, INPtables *tab)
     }
 #endif
 
+//#ifdef RELAN
+    if (strcmp (token, "relan") == 0)
+    {
+        JOB *relanJob ;
+        which = ft_find_analysis ("RELAN") ;
+        if (which == -1)
+        {
+            current->error = INPerrCat (current->error, INPmkTemp ("Reliability Analysis Unsupported\n")) ;
+            return (0) ; /* temporary */
+        }
+        err = ft_sim->newAnalysis (ft_curckt->ci_ckt, which, "relan", &relanJob, ft_curckt->ci_specTask) ;
+        if (err)
+        {
+            ft_sperror(err, "createRELAN") ;
+            return (0) ;
+        }
+
+        parm = INPgetValue (ckt, &line, IF_REAL, tab) ;
+        error = INPapName (ckt, which, relanJob, "relan_aging_step", parm) ;
+        if (error)
+        {
+            current->error = INPerrCat (current->error, INPerror (error)) ;
+        }
+
+        parm = INPgetValue (ckt, &line, IF_REAL, tab) ;
+        error = INPapName (ckt, which, relanJob, "relan_aging_stop", parm) ;
+        if (error)
+        {
+            current->error = INPerrCat (current->error, INPerror (error)) ;
+        }
+
+        if (*line)
+        {
+            if (*line == 'd')
+            {
+                goto next ;
+            } else if (*line == 'u') {
+                goto uic_relan ;
+            }
+            parm = INPgetValue (ckt, &line, IF_REAL, tab) ; /* AgingStart */
+            error = INPapName (ckt, which, relanJob, "relan_aging_start", parm) ;
+            if (error)
+            {
+                current->error = INPerrCat (current->error, INPerror (error)) ;
+            }
+
+        uic_relan:
+            if (*line == 'u')
+            {
+                INPgetTok (&line, &name, 1) ;
+                if (strcmp (name, "uic") == 0)
+                {
+                    ptemp.iValue = 1 ;
+                    error = INPapName (ckt, which, relanJob, "uic", &ptemp) ;
+                    if (error)
+                    {
+                        current->error = INPerrCat (current->error, INPerror (error)) ;
+                    }
+                }
+            }
+        }
+    }
+//#endif
+
 next:
     while (*line) { /* read the entire line */
 
diff --git a/src/frontend/spiceif.c b/src/frontend/spiceif.c
index 550559e53..74616d884 100644
--- a/src/frontend/spiceif.c
+++ b/src/frontend/spiceif.c
@@ -221,6 +221,11 @@ if_run(CKTcircuit *ckt, char *what, wordlist *args, INPtables *tab)
         || eq(what, "pss")
         /* SP */
 #endif
+
+//#ifdef RELAN
+        || eq (what, "relan")
+//#endif
+
         )
     {
         s = wl_flatten(args); /* va: tfree char's tmalloc'ed in wl_flatten */
@@ -336,6 +341,11 @@ if_run(CKTcircuit *ckt, char *what, wordlist *args, INPtables *tab)
         (eq(what, "pss")) ||
         /* SP */
 #endif
+
+//#ifdef RELAN
+        (eq (what, "relan")) ||
+//#endif
+
         (eq(what, "run")))
     {
         /*CDHW Run the analysis pointed to by ci_curTask CDHW*/
diff --git a/src/frontend/typesdef.c b/src/frontend/typesdef.c
index 71584e719..9263d230b 100644
--- a/src/frontend/typesdef.c
+++ b/src/frontend/typesdef.c
@@ -82,6 +82,11 @@ static struct plotab plotabs[NUMPLOTTYPES] = {
     { "harm", "harm" },
     { "spect", "spect" },
     { "pss", "periodic" },
+
+//#ifdef RELAN
+    { "relan", "reliability" },
+//#endif
+
 };
 
 
diff --git a/src/include/ngspice/Makefile.am b/src/include/ngspice/Makefile.am
index efd91ae32..6eb592a4d 100644
--- a/src/include/ngspice/Makefile.am
+++ b/src/include/ngspice/Makefile.am
@@ -125,6 +125,7 @@ pkginclude_HEADERS = \
 	FastNorm3.h	\
 	fftext.h	\
 	wallace.h	\
-	wstdio.h
+	wstdio.h	\
+	relandefs.h
 
 MAINTAINERCLEANFILES = Makefile.in
diff --git a/src/include/ngspice/cktdefs.h b/src/include/ngspice/cktdefs.h
index 96b52d2d3..afdc5d430 100644
--- a/src/include/ngspice/cktdefs.h
+++ b/src/include/ngspice/cktdefs.h
@@ -423,7 +423,13 @@ extern int DCpss(CKTcircuit *, int);
 /* SP */
 #endif
 
+//#ifdef RELAN
+extern int RELANaskQuest (CKTcircuit *, JOB *, int, IFvalue *) ;
+extern int RELANsetParm (CKTcircuit *, JOB *, int, IFvalue *) ;
+extern int RELANinit (CKTcircuit *, JOB *) ;
+extern int RELANanalysis (CKTcircuit *, int) ;
 extern int CKTreliability (CKTcircuit *, unsigned int) ;
+//#endif
 
 extern int NaskQuest(CKTcircuit *, JOB *, int, IFvalue *);
 extern int NsetParm(CKTcircuit *, JOB *, int, IFvalue *);
diff --git a/src/include/ngspice/relandefs.h b/src/include/ngspice/relandefs.h
new file mode 100644
index 000000000..22607b4ab
--- /dev/null
+++ b/src/include/ngspice/relandefs.h
@@ -0,0 +1,29 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#ifndef ngspice_RELANDEFS_H
+#define ngspice_RELANDEFS_H
+
+#include "ngspice/cktdefs.h"
+#include "ngspice/jobdefs.h"
+#include "ngspice/tskdefs.h"
+
+typedef struct {
+    int JOBtype ;
+    JOB *JOBnextJob ;
+    char *JOBname ;
+    double RELANfinalTime ;
+    double RELANstep ;
+    double RELANmaxStep ;
+    double RELANinitTime ;
+    long RELANmode ;
+    runDesc *RELANplot ;
+} RELANan ;
+
+#define RELAN_TSTART 1
+#define RELAN_TSTOP 2
+#define RELAN_TSTEP 3
+#define RELAN_TMAX 4
+#define RELAN_UIC 5
+#endif
diff --git a/src/spicelib/analysis/Makefile.am b/src/spicelib/analysis/Makefile.am
index 2cdceb11b..6d10c779f 100644
--- a/src/spicelib/analysis/Makefile.am
+++ b/src/spicelib/analysis/Makefile.am
@@ -96,7 +96,11 @@ libckt_la_SOURCES = \
 		traninit.c	\
 		transetp.c	\
 		cluster.c	\
-		cktreliability.c
+		cktreliability.c	\
+		relananalysis.c	\
+		relanaskq.c	\
+		relaninit.c	\
+		relansetp.c
 
 
 if PSS_WANTED
diff --git a/src/spicelib/analysis/analysis.c b/src/spicelib/analysis/analysis.c
index ecb9279a6..fb7f3dfcd 100644
--- a/src/spicelib/analysis/analysis.c
+++ b/src/spicelib/analysis/analysis.c
@@ -5,27 +5,29 @@
 
 #include "analysis.h"
 
-extern SPICEanalysis OPTinfo;
-extern SPICEanalysis ACinfo;
-extern SPICEanalysis DCTinfo;
-extern SPICEanalysis DCOinfo;
-extern SPICEanalysis TRANinfo;
-extern SPICEanalysis PZinfo;
-extern SPICEanalysis TFinfo;
-extern SPICEanalysis DISTOinfo;
-extern SPICEanalysis NOISEinfo;
-extern SPICEanalysis SENSinfo;
+extern SPICEanalysis OPTinfo ;
+extern SPICEanalysis RELANinfo ;
+extern SPICEanalysis ACinfo ;
+extern SPICEanalysis DCTinfo ;
+extern SPICEanalysis DCOinfo ;
+extern SPICEanalysis TRANinfo ;
+extern SPICEanalysis PZinfo ;
+extern SPICEanalysis TFinfo ;
+extern SPICEanalysis DISTOinfo ;
+extern SPICEanalysis NOISEinfo ;
+extern SPICEanalysis SENSinfo ;
 
 #ifdef WITH_PSS
-extern SPICEanalysis PSSinfo;
+extern SPICEanalysis PSSinfo ;
 #endif
 
 #ifdef WANT_SENSE2
-extern SPICEanalysis SEN2info;
+extern SPICEanalysis SEN2info ;
 #endif
 
-SPICEanalysis *analInfo[] = {
+SPICEanalysis *analInfo [] = {
     &OPTinfo,
+    &RELANinfo,
     &ACinfo,
     &DCTinfo,
     &DCOinfo,
@@ -41,25 +43,24 @@ SPICEanalysis *analInfo[] = {
 #ifdef WANT_SENSE2
     &SEN2info,
 #endif
-};
+} ;
 
-
-char *spice_analysis_get_name(int index)
+char *spice_analysis_get_name (int index)
 {
-    return analInfo[index]->if_analysis.name;
+    return analInfo [index]->if_analysis.name ;
 }
 
-char *spice_analysis_get_description(int index)
+char *spice_analysis_get_description (int index)
 {
-    return analInfo[index]->if_analysis.description;
+    return analInfo [index]->if_analysis.description ;
 }
 
-int spice_num_analysis(void)
+int spice_num_analysis (void)
 {
-    return NUMELEMS(analInfo);
+    return NUMELEMS (analInfo) ;
 }
 
-SPICEanalysis **spice_analysis_ptr(void)
+SPICEanalysis **spice_analysis_ptr (void)
 {
-    return analInfo;
+    return analInfo ;
 }
diff --git a/src/spicelib/analysis/cktdojob.c b/src/spicelib/analysis/cktdojob.c
index 93cbd3b39..72e755b02 100644
--- a/src/spicelib/analysis/cktdojob.c
+++ b/src/spicelib/analysis/cktdojob.c
@@ -224,7 +224,17 @@ CKTdoJob(CKTcircuit *ckt, int reset, TSKtask *task)
 			/* txl, cpl addition */
 			if (error == 1111) break;
 
-                    CKTreliability (ckt, 1) ;
+                    if (i == 1)
+                    {
+                        /* In case of Reliability Analysis, perform the final CKTreliability */
+                        printf ("\n\nFinal Aging...\n") ;
+                        CKTreliability (ckt, 1) ;
+                        error = CKTtemp (ckt) ;
+                        if (error)
+                        {
+                            return (error) ;
+                        }
+                    }
 
 		}
 		if (error)
diff --git a/src/spicelib/analysis/dctran.c b/src/spicelib/analysis/dctran.c
index 3d0b1d263..f43189898 100644
--- a/src/spicelib/analysis/dctran.c
+++ b/src/spicelib/analysis/dctran.c
@@ -828,17 +828,6 @@ resume:
 #endif
 
         } else {
-
-            printf ("CKTtime: %-.9g\n", ckt->CKTtime) ;
-            CKTreliability (ckt, 0) ;
-/** Nel RHSold, ogni device deve accedere ai propri valori per vedere se esso stesso è acceso o spento.
-    Nel caso del BSIM4, vale la regola Vgs > Vth, dove Vgs = ckt->CKTrhsOld [here->BSIM4...] - ckt->CKTrhsOld [here->BSIM4...] e Vth = here->BSIM4vth .
-    In caso il transistor sia acceso, si alza un flag, privato del device, che indica che è acceso. Se è spento, lo stesso flag sarà basso.
-    Il tempo corrente CKTtime deve essere memorizzato insieme, in modo tale da poter poi calcolare il delta di tempo necessario al modello.
-    QUI, deve essere controllato che all'istante precedente il device sia acceso (o spento). Se si manifesta un cambio, allora la fase di stress (o di recovery)
-      è finita e bisogna calcolare il delta_vth attraverso il modello.
-*/
-
             if (firsttime) {
 #ifdef WANT_SENSE2
                 if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
diff --git a/src/spicelib/analysis/relananalysis.c b/src/spicelib/analysis/relananalysis.c
new file mode 100644
index 000000000..bef3d9ce4
--- /dev/null
+++ b/src/spicelib/analysis/relananalysis.c
@@ -0,0 +1,994 @@
+/**********
+Copyright 1990 Regents of the University of California.  All rights reserved.
+Author: 1985 Thomas L. Quarles
+Modified: 2000  AlansFixes
+**********/
+
+/* subroutine to do DC TRANSIENT analysis
+        --- ONLY, unlike spice2 routine with the same name! */
+
+#include "ngspice/ngspice.h"
+#include "ngspice/cktdefs.h"
+#include "cktaccept.h"
+#include "ngspice/relandefs.h"
+#include "ngspice/sperror.h"
+#include "ngspice/fteext.h"
+#include "ngspice/missing_math.h"
+
+/* for setting breakpoints required by dbs data base */
+extern struct dbcomm *dbs;
+#include "ngspice/ftedebug.h"
+
+#ifdef XSPICE
+/* gtri - add - wbk - Add headers */
+#include "ngspice/miftypes.h"
+
+#include "ngspice/evt.h"
+#include "ngspice/mif.h"
+#include "ngspice/evtproto.h"
+#include "ngspice/ipctiein.h"
+/* gtri - end - wbk - Add headers */
+#endif
+
+#ifdef CLUSTER
+#include "ngspice/cluster.h"
+#endif
+
+#ifdef SHARED_MODULE
+extern int add_bkpt(void);
+extern int sharedsync(double*, double*, double, double, double, int, int*, int);
+extern int ng_ident;      /* for debugging */
+static double del_before; /* for debugging */
+#endif
+
+#define INIT_STATS() \
+do { \
+    startTime = SPfrontEnd->IFseconds();        \
+    startIters = ckt->CKTstat->STATnumIter;     \
+    startdTime = ckt->CKTstat->STATdecompTime;  \
+    startsTime = ckt->CKTstat->STATsolveTime;   \
+    startlTime = ckt->CKTstat->STATloadTime;    \
+    startkTime = ckt->CKTstat->STATsyncTime;    \
+} while(0)
+
+#define UPDATE_STATS(analysis) \
+do { \
+    ckt->CKTcurrentAnalysis = analysis; \
+    ckt->CKTstat->STATtranTime += SPfrontEnd->IFseconds() - startTime; \
+    ckt->CKTstat->STATtranIter += ckt->CKTstat->STATnumIter - startIters; \
+    ckt->CKTstat->STATtranDecompTime += ckt->CKTstat->STATdecompTime - startdTime; \
+    ckt->CKTstat->STATtranSolveTime += ckt->CKTstat->STATsolveTime - startsTime; \
+    ckt->CKTstat->STATtranLoadTime += ckt->CKTstat->STATloadTime - startlTime; \
+    ckt->CKTstat->STATtranSyncTime += ckt->CKTstat->STATsyncTime - startkTime; \
+} while(0)
+
+
+int
+RELANanalysis (CKTcircuit *ckt,
+       int restart)   /* forced restart flag */
+{
+    RELANan *job = (RELANan *) ckt->CKTcurJob;
+
+    int i;
+    double olddelta;
+    double delta;
+    double newdelta;
+    double *temp;
+    double startdTime;
+    double startsTime;
+    double startlTime;
+    double startkTime;
+    double startTime;
+    int startIters;
+    int converged;
+    int firsttime;
+    int error;
+#ifdef WANT_SENSE2
+    int save, save2, size;
+    long save1;
+#endif
+    int save_order;
+    long save_mode;
+    IFuid timeUid;
+    IFuid *nameList;
+    int numNames;
+    double maxstepsize = 0.0;
+
+    int ltra_num;
+    CKTnode *node;
+#ifdef XSPICE
+/* gtri - add - wbk - 12/19/90 - Add IPC stuff */
+    Ipc_Boolean_t  ipc_firsttime = IPC_TRUE;
+    Ipc_Boolean_t  ipc_secondtime = IPC_FALSE;
+    Ipc_Boolean_t  ipc_delta_cut = IPC_FALSE;
+    double         ipc_last_time = 0.0;
+    double         ipc_last_delta = 0.0;
+/* gtri - end - wbk - 12/19/90 - Add IPC stuff */
+#endif
+#if defined CLUSTER || defined SHARED_MODULE
+    int redostep;
+#endif
+    if(restart || ckt->CKTtime == 0) {
+        delta=MIN(ckt->CKTfinalTime/100,ckt->CKTstep)/10;
+
+#ifdef STEPDEBUG
+        printf("delta = %g    finalTime/100: %g    CKTstep: %g\n",delta,ckt->CKTfinalTime/100,ckt->CKTstep);
+#endif
+        /* begin LTRA code addition */
+        if (ckt->CKTtimePoints != NULL)
+            FREE(ckt->CKTtimePoints);
+
+        if (ckt->CKTstep >= ckt->CKTmaxStep)
+            maxstepsize = ckt->CKTstep;
+        else
+            maxstepsize = ckt->CKTmaxStep;
+
+        ckt->CKTsizeIncr = 10;
+        ckt->CKTtimeIndex = -1; /* before the DC soln has been stored */
+        ckt->CKTtimeListSize = (int) ceil( ckt->CKTfinalTime / maxstepsize );
+        ltra_num = CKTtypelook("LTRA");
+        if (ltra_num >= 0 && ckt->CKThead[ltra_num] != NULL)
+            ckt->CKTtimePoints = TMALLOC(double, ckt->CKTtimeListSize);
+        /* end LTRA code addition */
+
+        if(ckt->CKTbreaks) FREE(ckt->CKTbreaks);
+        ckt->CKTbreaks = TMALLOC(double, 2);
+        if(ckt->CKTbreaks == NULL) return(E_NOMEM);
+        ckt->CKTbreaks[0] = 0;
+        ckt->CKTbreaks[1] = ckt->CKTfinalTime;
+        ckt->CKTbreakSize = 2;
+
+#ifdef SHARED_MODULE
+        add_bkpt();
+#endif
+
+#ifdef XSPICE
+/* gtri - begin - wbk - 12/19/90 - Modify setting of CKTminBreak */
+        /* Set to 10 times delmin for ATESSE 1 compatibity */
+        if(ckt->CKTminBreak==0) ckt->CKTminBreak = 10.0 * ckt->CKTdelmin;
+/* gtri - end - wbk - 12/19/90 - Modify setting of CKTminBreak */
+#else
+        if(ckt->CKTminBreak==0) ckt->CKTminBreak=ckt->CKTmaxStep*5e-5;
+#endif
+
+#ifdef XSPICE
+/* gtri - add - wbk - 12/19/90 - Add IPC stuff and set anal_init and anal_type */
+        /* Tell the beginPlot routine what mode we're in */
+        g_ipc.anal_type = IPC_ANAL_TRAN;
+
+        /* Tell the code models what mode we're in */
+        g_mif_info.circuit.anal_type = MIF_DC;
+
+        g_mif_info.circuit.anal_init = MIF_TRUE;
+/* gtri - end - wbk */
+#endif
+        error = CKTnames(ckt,&numNames,&nameList);
+        if(error) return(error);
+        SPfrontEnd->IFnewUid (ckt, &timeUid, NULL, "time", UID_OTHER, NULL);
+        error = SPfrontEnd->OUTpBeginPlot (ckt, ckt->CKTcurJob,
+                                           ckt->CKTcurJob->JOBname,
+                                           timeUid, IF_REAL,
+                                           numNames, nameList, IF_REAL,
+                                           &(job->RELANplot));
+        tfree(nameList);
+        if(error) return(error);
+
+        /* initialize CKTsoaCheck `warn' counters */
+        if (ckt->CKTsoaCheck)
+            error = CKTsoaInit();
+
+        ckt->CKTtime = 0;
+        ckt->CKTdelta = 0;
+        ckt->CKTbreak = 1;
+        firsttime = 1;
+        save_mode = (ckt->CKTmode&MODEUIC) | MODETRANOP | MODEINITJCT;
+        save_order = ckt->CKTorder;
+
+/* Add breakpoints here which have been requested by the user setting the
+   stop command as 'stop when time = xx'.
+   Get data from the global dbs data base.
+*/
+        if (dbs) {
+            struct dbcomm *d;
+            for (d = dbs; d; d = d->db_next)
+                if ((d->db_type == DB_STOPWHEN) && cieq(d->db_nodename1,"time")
+                    && (d->db_value2 > 0)) {
+                    CKTsetBreak(ckt, d->db_value2);
+                    if (ft_ngdebug)
+                        printf("breakpoint set to time = %g\n", d->db_value2);
+                }
+        }
+
+#ifdef XSPICE
+/* gtri - begin - wbk - set a breakpoint at end of supply ramping time */
+        /* must do this after CKTtime set to 0 above */
+        if(ckt->enh->ramp.ramptime > 0.0)
+            CKTsetBreak(ckt, ckt->enh->ramp.ramptime);
+/* gtri - end - wbk - set a breakpoint at end of supply ramping time */
+
+/* gtri - begin - wbk - Call EVTop if event-driven instances exist */
+        if(ckt->evt->counts.num_insts != 0) {
+            /* use new DCOP algorithm */
+            converged = EVTop(ckt,
+                        (ckt->CKTmode & MODEUIC) | MODETRANOP | MODEINITJCT,
+                        (ckt->CKTmode & MODEUIC) | MODETRANOP | MODEINITFLOAT,
+                        ckt->CKTdcMaxIter,
+                        MIF_TRUE);
+            EVTdump(ckt, IPC_ANAL_DCOP, 0.0);
+
+            EVTop_save(ckt, MIF_FALSE, 0.0);
+
+/* gtri - end - wbk - Call EVTop if event-driven instances exist */
+        } else
+#endif
+            converged = CKTop(ckt,
+                (ckt->CKTmode & MODEUIC) | MODETRANOP | MODEINITJCT,
+                (ckt->CKTmode & MODEUIC) | MODETRANOP | MODEINITFLOAT,
+                ckt->CKTdcMaxIter);
+
+        if(converged != 0) {
+            fprintf(stdout,"\nTransient solution failed -\n");
+            CKTncDump(ckt);
+/*          CKTnode *node;
+            double new, old, tol;
+            int i=1;
+
+            fprintf(stdout,"\nTransient solution failed -\n\n");
+            fprintf(stdout,"Last Node Voltages\n");
+            fprintf(stdout,"------------------\n\n");
+            fprintf(stdout,"%-30s %20s %20s\n", "Node", "Last Voltage",
+                                                               "Previous Iter");
+            fprintf(stdout,"%-30s %20s %20s\n", "----", "------------",
+                                                               "-------------");
+            for(node=ckt->CKTnodes->next;node;node=node->next) {
+                if (strstr(node->name, "#branch") || !strstr(node->name, "#")) {
+                    new =  ckt->CKTrhsOld [i] ;
+                    old =  ckt->CKTrhs [i] ;
+                    fprintf(stdout,"%-30s %20g %20g", node->name, new, old);
+                    if(node->type == SP_VOLTAGE) {
+                        tol =  ckt->CKTreltol * (MAX(fabs(old),fabs(new))) +
+                                ckt->CKTvoltTol;
+                    } else {
+                        tol =  ckt->CKTreltol * (MAX(fabs(old),fabs(new))) +
+                                ckt->CKTabstol;
+                    }
+                    if (fabs(new-old) >tol ) {
+                        fprintf(stdout," *");
+                    }
+                    fprintf(stdout,"\n");
+                }
+                i++;
+            } */
+            fprintf(stdout,"\n");
+            fflush(stdout);
+        } else if (!ft_noacctprint && !ft_noinitprint) {
+            fprintf(stdout,"\nInitial Transient Solution\n");
+            fprintf(stdout,"--------------------------\n\n");
+            fprintf(stdout,"%-30s %15s\n", "Node", "Voltage");
+            fprintf(stdout,"%-30s %15s\n", "----", "-------");
+            for(node=ckt->CKTnodes->next;node;node=node->next) {
+                if (strstr(node->name, "#branch") || !strstr(node->name, "#"))
+                    fprintf(stdout,"%-30s %15g\n", node->name,
+                                              ckt->CKTrhsOld[node->number]);
+            }
+            fprintf(stdout,"\n");
+            fflush(stdout);
+        }
+
+        if(converged != 0) return(converged);
+#ifdef XSPICE
+/* gtri - add - wbk - 12/19/90 - Add IPC stuff */
+
+        /* Send the operating point results for Mspice compatibility */
+        if(g_ipc.enabled) {
+            ipc_send_dcop_prefix();
+            CKTdump(ckt, 0.0, job->RELANplot);
+            ipc_send_dcop_suffix();
+        }
+
+/* gtri - end - wbk */
+
+/* gtri - add - wbk - 12/19/90 - set anal_init and anal_type */
+
+        g_mif_info.circuit.anal_init = MIF_TRUE;
+
+        /* Tell the code models what mode we're in */
+        g_mif_info.circuit.anal_type = MIF_TRAN;
+
+/* gtri - end - wbk */
+
+/* gtri - begin - wbk - Add Breakpoint stuff */
+
+        /* Initialize the temporary breakpoint variables to infinity */
+        g_mif_info.breakpoint.current = 1.0e30;
+        g_mif_info.breakpoint.last    = 1.0e30;
+
+/* gtri - end - wbk - Add Breakpoint stuff */
+#endif
+        ckt->CKTstat->STATtimePts ++;
+        ckt->CKTorder = 1;
+        for(i=0;i<7;i++) {
+            ckt->CKTdeltaOld[i]=ckt->CKTmaxStep;
+        }
+        ckt->CKTdelta = delta;
+#ifdef STEPDEBUG
+        (void)printf("delta initialized to %g\n",ckt->CKTdelta);
+#endif
+        ckt->CKTsaveDelta = ckt->CKTfinalTime/50;
+
+#ifdef WANT_SENSE2
+        if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
+#ifdef SENSDEBUG
+            printf("\nTransient Sensitivity Results\n\n");
+            CKTsenPrint(ckt);
+#endif /* SENSDEBUG */
+            save = ckt->CKTsenInfo->SENmode;
+            ckt->CKTsenInfo->SENmode = TRANSEN;
+            save1 = ckt->CKTmode;
+            save2 = ckt->CKTorder;
+            ckt->CKTmode = save_mode;
+            ckt->CKTorder = save_order;
+            error = CKTsenDCtran(ckt);
+            if (error)
+                return(error);
+
+            ckt->CKTmode = save1;
+            ckt->CKTorder = save2;
+        }
+#endif
+
+        ckt->CKTmode = (ckt->CKTmode&MODEUIC) | MODETRAN | MODEINITTRAN;
+        /* modeinittran set here */
+        ckt->CKTag[0]=ckt->CKTag[1]=0;
+        bcopy(ckt->CKTstate0, ckt->CKTstate1,
+              (size_t) ckt->CKTnumStates * sizeof(double));
+
+#ifdef WANT_SENSE2
+        if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
+            size = SMPmatSize(ckt->CKTmatrix);
+            for(i = 1; i<=size ; i++)
+                ckt->CKTrhsOp[i] = ckt->CKTrhsOld[i];
+        }
+#endif
+
+        INIT_STATS();
+#ifdef CLUSTER
+        CLUsetup(ckt);
+#endif
+    } else {
+        /* saj As traninit resets CKTmode */
+        ckt->CKTmode = (ckt->CKTmode&MODEUIC) | MODETRAN | MODEINITPRED;
+        /* saj */
+        INIT_STATS();
+        if(ckt->CKTminBreak==0) ckt->CKTminBreak=ckt->CKTmaxStep*5e-5;
+        firsttime=0;
+        /* To get rawfile working saj*/
+        error = SPfrontEnd->OUTpBeginPlot (NULL, NULL,
+                                           NULL,
+                                           NULL, 0,
+                                           666, NULL, 666,
+                                           &(job->RELANplot));
+        if(error) {
+            fprintf(stderr, "Couldn't relink rawfile\n");
+            return error;
+        }
+        /* end saj*/
+        goto resume;
+    }
+
+/* 650 */
+    nextTime:
+
+    /* begin LTRA code addition */
+    if (ckt->CKTtimePoints) {
+    ckt->CKTtimeIndex++;
+        if (ckt->CKTtimeIndex >= ckt->CKTtimeListSize) {
+            /* need more space */
+            int need;
+            need = (int) ceil( (ckt->CKTfinalTime - ckt->CKTtime) / maxstepsize );
+            if (need < ckt->CKTsizeIncr)
+                need = ckt->CKTsizeIncr;
+            ckt->CKTtimeListSize += need;
+            ckt->CKTtimePoints = TREALLOC(double, ckt->CKTtimePoints, ckt->CKTtimeListSize);
+            ckt->CKTsizeIncr = (int) ceil(1.4 * ckt->CKTsizeIncr);
+        }
+        ckt->CKTtimePoints[ckt->CKTtimeIndex] = ckt->CKTtime;
+    }
+    /* end LTRA code addition */
+
+    error = CKTaccept(ckt);
+    /* check if current breakpoint is outdated; if so, clear */
+    if (ckt->CKTtime > ckt->CKTbreaks[0]) CKTclrBreak(ckt);
+
+    if (ckt->CKTsoaCheck)
+        error = CKTsoaCheck(ckt);
+
+    /*
+ * Breakpoint handling scheme:
+ * When a timepoint t is accepted (by CKTaccept), clear all previous
+ * breakpoints, because they will never be needed again.
+ *
+ * t may itself be a breakpoint, or indistinguishably close. DON'T
+ * clear t itself; recognise it as a breakpoint and act accordingly
+ *
+ * if t is not a breakpoint, limit the timestep so that the next
+ * breakpoint is not crossed
+ */
+
+#ifdef STEPDEBUG
+    printf("Delta %g accepted at time %g (finaltime: %g)\n",ckt->CKTdelta,ckt->CKTtime,ckt->CKTfinalTime);
+    fflush(stdout);
+#endif /* STEPDEBUG */
+    ckt->CKTstat->STATaccepted ++;
+    ckt->CKTbreak = 0;
+    /* XXX Error will cause single process to bail. */
+    if(error)  {
+        UPDATE_STATS(DOING_TRAN);
+        return(error);
+    }
+#ifdef XSPICE
+/* gtri - modify - wbk - 12/19/90 - Send IPC stuff */
+
+    if(g_ipc.enabled) {
+
+        /* Send event-driven results */
+        EVTdump(ckt, IPC_ANAL_TRAN, 0.0);
+
+        /* Then follow with analog results... */
+
+        /* Test to see if delta was cut by a breakpoint, */
+        /* a non-convergence, or a too large truncation error */
+        if(ipc_firsttime)
+            ipc_delta_cut = IPC_FALSE;
+        else if(ckt->CKTtime < (ipc_last_time + (0.999 * ipc_last_delta)))
+            ipc_delta_cut = IPC_TRUE;
+        else
+            ipc_delta_cut = IPC_FALSE;
+
+        /* Record the data required to check for delta cuts */
+        ipc_last_time = ckt->CKTtime;
+        ipc_last_delta = MIN(ckt->CKTdelta, ckt->CKTmaxStep);
+
+        /* Send results data if time since last dump is greater */
+        /* than 'mintime', or if first or second timepoints, */
+        /* or if delta was cut */
+        if( (ckt->CKTtime >= (g_ipc.mintime + g_ipc.last_time)) ||
+            ipc_firsttime || ipc_secondtime || ipc_delta_cut ) {
+
+            ipc_send_data_prefix(ckt->CKTtime);
+            CKTdump(ckt, ckt->CKTtime, job->RELANplot);
+            ipc_send_data_suffix();
+
+            if(ipc_firsttime) {
+                ipc_firsttime = IPC_FALSE;
+                ipc_secondtime = IPC_TRUE;
+            } else if(ipc_secondtime) {
+                ipc_secondtime = IPC_FALSE;
+            }
+
+            g_ipc.last_time = ckt->CKTtime;
+        }
+    } else
+/* gtri - modify - wbk - 12/19/90 - Send IPC stuff */
+#endif
+#ifdef CLUSTER
+        CLUoutput(ckt);
+#endif
+        if(ckt->CKTtime >= ckt->CKTinitTime)
+            CKTdump(ckt, ckt->CKTtime, job->RELANplot);
+#ifdef XSPICE
+/* gtri - begin - wbk - Update event queues/data for accepted timepoint */
+    /* Note: this must be done AFTER sending results to SI so it can't */
+    /* go next to CKTaccept() above */
+    if(ckt->evt->counts.num_insts > 0)
+        EVTaccept(ckt, ckt->CKTtime);
+/* gtri - end - wbk - Update event queues/data for accepted timepoint */
+#endif
+    ckt->CKTstat->STAToldIter = ckt->CKTstat->STATnumIter;
+    if(check_autostop("tran") ||
+      fabs(ckt->CKTtime - ckt->CKTfinalTime) < ckt->CKTminBreak ||
+      AlmostEqualUlps( ckt->CKTtime, ckt->CKTfinalTime, 100 ) ) {
+#ifdef STEPDEBUG
+        printf(" done:  time is %g, final time is %g, and tol is %g\n",
+        ckt->CKTtime, ckt->CKTfinalTime, ckt->CKTminBreak);
+#endif
+        SPfrontEnd->OUTendPlot (job->RELANplot);
+        job->RELANplot = NULL;
+        UPDATE_STATS(0);
+#ifdef WANT_SENSE2
+        if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
+            ckt->CKTsenInfo->SENmode = save;
+        }
+#endif
+        return(OK);
+    }
+    if(SPfrontEnd->IFpauseTest()) {
+        /* user requested pause... */
+        UPDATE_STATS(DOING_TRAN);
+        return(E_PAUSE);
+    }
+resume:
+#ifdef STEPDEBUG
+    if( (ckt->CKTdelta <= ckt->CKTfinalTime/50) &&
+        (ckt->CKTdelta <= ckt->CKTmaxStep)) {
+        ;
+    } else {
+        if(ckt->CKTfinalTime/50<ckt->CKTmaxStep) {
+            (void)printf("limited by Tstop/50\n");
+        } else {
+            (void)printf("limited by Tmax == %g\n",ckt->CKTmaxStep);
+        }
+    }
+#endif
+#ifdef HAS_PROGREP
+    if (ckt->CKTtime == 0.)
+        SetAnalyse( "tran init", 0);
+    else
+        SetAnalyse( "tran", (int)((ckt->CKTtime * 1000.) / ckt->CKTfinalTime + 0.5));
+#endif
+    ckt->CKTdelta =
+            MIN(ckt->CKTdelta,ckt->CKTmaxStep);
+#ifdef XSPICE
+/* gtri - begin - wbk - Cut integration order if first timepoint after breakpoint */
+    /* if(ckt->CKTtime == g_mif_info.breakpoint.last) */
+    if ( AlmostEqualUlps( ckt->CKTtime, g_mif_info.breakpoint.last, 100 ) )
+        ckt->CKTorder = 1;
+/* gtri - end   - wbk - Cut integration order if first timepoint after breakpoint */
+
+#endif
+
+  /* are we at a breakpoint, or indistinguishably close? */
+    /* if ((ckt->CKTtime == ckt->CKTbreaks[0]) || (ckt->CKTbreaks[0] - */
+    if ( AlmostEqualUlps( ckt->CKTtime, ckt->CKTbreaks[0], 100 ) ||
+         ckt->CKTbreaks[0] - ckt->CKTtime <= ckt->CKTdelmin) {
+        /* first timepoint after a breakpoint - cut integration order */
+        /* and limit timestep to .1 times minimum of time to next breakpoint,
+         * and previous timestep
+         */
+        ckt->CKTorder = 1;
+#ifdef STEPDEBUG
+        if( (ckt->CKTdelta > .1*ckt->CKTsaveDelta) ||
+            (ckt->CKTdelta > .1*(ckt->CKTbreaks[1] - ckt->CKTbreaks[0])) ) {
+            if(ckt->CKTsaveDelta < (ckt->CKTbreaks[1] - ckt->CKTbreaks[0]))  {
+                (void)printf("limited by pre-breakpoint delta (saveDelta: %g, nxt_breakpt: %g, curr_breakpt: %g\n",
+                  ckt->CKTsaveDelta, ckt->CKTbreaks[1], ckt->CKTbreaks[0]);
+            } else {
+                (void)printf("limited by next breakpoint\n");
+            }
+        }
+#endif
+
+        ckt->CKTdelta = MIN(ckt->CKTdelta, .1 * MIN(ckt->CKTsaveDelta,
+            ckt->CKTbreaks[1] - ckt->CKTbreaks[0]));
+
+        if(firsttime) {
+            /* set a breakpoint to reduce ringing of current in devices */
+            if (ckt->CKTmode & MODEUIC)
+                CKTsetBreak(ckt, ckt->CKTstep);
+
+            ckt->CKTdelta /= 10;
+#ifdef STEPDEBUG
+            (void)printf("delta cut for initial timepoint\n");
+#endif
+        }
+
+#ifndef XSPICE
+        /* don't want to get below delmin for no reason */
+        ckt->CKTdelta = MAX(ckt->CKTdelta, ckt->CKTdelmin*2.0);
+#endif
+
+    }
+
+#ifndef XSPICE
+    else if(ckt->CKTtime + ckt->CKTdelta >= ckt->CKTbreaks[0]) {
+        ckt->CKTsaveDelta = ckt->CKTdelta;
+        ckt->CKTdelta = ckt->CKTbreaks[0] - ckt->CKTtime;
+#ifdef STEPDEBUG
+        (void)printf("delta cut to %g to hit breakpoint\n",ckt->CKTdelta);
+        fflush(stdout);
+#endif
+        ckt->CKTbreak = 1; /* why? the current pt. is not a bkpt. */
+    }
+#endif /* !XSPICE */
+
+
+#ifdef XSPICE
+/* gtri - begin - wbk - Add Breakpoint stuff */
+
+    if(ckt->CKTtime + ckt->CKTdelta >= g_mif_info.breakpoint.current) {
+        /* If next time > temporary breakpoint, force it to the breakpoint */
+        /* And mark that timestep was set by temporary breakpoint */
+        ckt->CKTsaveDelta = ckt->CKTdelta;
+        ckt->CKTdelta = g_mif_info.breakpoint.current - ckt->CKTtime;
+        g_mif_info.breakpoint.last = ckt->CKTtime + ckt->CKTdelta;
+    } else {
+        /* Else, mark that timestep was not set by temporary breakpoint */
+        g_mif_info.breakpoint.last = 1.0e30;
+    }
+
+/* gtri - end - wbk - Add Breakpoint stuff */
+
+/* gtri - begin - wbk - Modify Breakpoint stuff */
+    /* Throw out any permanent breakpoint times <= current time */
+    for (;;) {
+#ifdef STEPDEBUG
+        printf("    brk_pt: %g    ckt_time: %g    ckt_min_break: %g\n",ckt->CKTbreaks[0], ckt->CKTtime, ckt->CKTminBreak);
+#endif
+        if(AlmostEqualUlps(ckt->CKTbreaks[0], ckt->CKTtime, 100) ||
+           ckt->CKTbreaks[0] <= ckt->CKTtime + ckt->CKTminBreak) {
+#ifdef STEPDEBUG
+            printf("throwing out permanent breakpoint times <= current time (brk pt: %g)\n",ckt->CKTbreaks[0]);
+            printf("    ckt_time: %g    ckt_min_break: %g\n",ckt->CKTtime, ckt->CKTminBreak);
+#endif
+            CKTclrBreak(ckt);
+        } else {
+            break;
+        }
+    }
+    /* Force the breakpoint if appropriate */
+    if(ckt->CKTtime + ckt->CKTdelta > ckt->CKTbreaks[0]) {
+        ckt->CKTbreak = 1;
+        ckt->CKTsaveDelta = ckt->CKTdelta;
+        ckt->CKTdelta = ckt->CKTbreaks[0] - ckt->CKTtime;
+    }
+
+/* gtri - end - wbk - Modify Breakpoint stuff */
+
+#ifdef SHARED_MODULE
+        /* Either directly go to next time step, or modify ckt->CKTdelta depending on
+           synchronization requirements. sharedsync() returns 0. */
+    sharedsync(&ckt->CKTtime, &ckt->CKTdelta, 0, ckt->CKTfinalTime,
+        ckt->CKTdelmin, 0, &ckt->CKTstat->STATrejected, 0);
+#endif
+
+/* gtri - begin - wbk - Do event solution */
+
+    if(ckt->evt->counts.num_insts > 0) {
+
+        /* if time = 0 and op_alternate was specified as false during */
+        /* dcop analysis, call any changed instances to let them */
+        /* post their outputs with their associated delays */
+        if((ckt->CKTtime == 0.0) && (! ckt->evt->options.op_alternate))
+            EVTiter(ckt);
+
+        /* while there are events on the queue with event time <= next */
+        /* projected analog time, process them */
+        while((g_mif_info.circuit.evt_step = EVTnext_time(ckt))
+               <= (ckt->CKTtime + ckt->CKTdelta)) {
+
+            /* Initialize temp analog bkpt to infinity */
+            g_mif_info.breakpoint.current = 1e30;
+
+            /* Pull items off queue and process them */
+            EVTdequeue(ckt, g_mif_info.circuit.evt_step);
+            EVTiter(ckt);
+
+            /* If any instances have forced an earlier */
+            /* next analog time, cut the delta */
+            if(ckt->CKTbreaks[0] < g_mif_info.breakpoint.current)
+                if(ckt->CKTbreaks[0] > ckt->CKTtime + ckt->CKTminBreak)
+                    g_mif_info.breakpoint.current = ckt->CKTbreaks[0];
+            if(g_mif_info.breakpoint.current < ckt->CKTtime + ckt->CKTdelta) {
+                /* Breakpoint must be > last accepted timepoint */
+                /* and >= current event time */
+                if(g_mif_info.breakpoint.current >  ckt->CKTtime + ckt->CKTminBreak  &&
+                   g_mif_info.breakpoint.current >= g_mif_info.circuit.evt_step) {
+                    ckt->CKTsaveDelta = ckt->CKTdelta;
+                    ckt->CKTdelta = g_mif_info.breakpoint.current - ckt->CKTtime;
+                    g_mif_info.breakpoint.last = ckt->CKTtime + ckt->CKTdelta;
+                }
+            }
+
+        } /* end while next event time <= next analog time */
+    } /* end if there are event instances */
+
+/* gtri - end - wbk - Do event solution */
+#else
+
+#ifdef CLUSTER
+    if(!CLUsync(ckt->CKTtime,&ckt->CKTdelta,0)) {
+      printf("Sync error!\n");
+      exit(0);
+    }
+#endif /* CLUSTER */
+
+#ifdef SHARED_MODULE
+    /* Either directly go to next time step, or modify ckt->CKTdelta depending on
+       synchronization requirements. sharedsync() returns 0.
+    */
+    sharedsync(&ckt->CKTtime, &ckt->CKTdelta, 0, ckt->CKTfinalTime,
+        ckt->CKTdelmin, 0, &ckt->CKTstat->STATrejected, 0);
+#endif
+
+#endif
+    for(i=5; i>=0; i--)
+        ckt->CKTdeltaOld[i+1] = ckt->CKTdeltaOld[i];
+    ckt->CKTdeltaOld[0] = ckt->CKTdelta;
+
+    temp = ckt->CKTstates[ckt->CKTmaxOrder+1];
+    for(i=ckt->CKTmaxOrder;i>=0;i--) {
+        ckt->CKTstates[i+1] = ckt->CKTstates[i];
+    }
+    ckt->CKTstates[0] = temp;
+
+/* 600 */
+    for (;;) {
+#if defined CLUSTER || defined SHARED_MODULE
+        redostep = 1;
+#endif
+#ifdef XSPICE
+/* gtri - add - wbk - 4/17/91 - Fix Berkeley bug */
+/* This is needed here to allow CAPask to output currents */
+/* during Transient analysis.  A grep for CKTcurrentAnalysis */
+/* indicates that it should not hurt anything else ... */
+
+        ckt->CKTcurrentAnalysis = DOING_TRAN;
+
+/* gtri - end - wbk - 4/17/91 - Fix Berkeley bug */
+#endif
+        olddelta=ckt->CKTdelta;
+        /* time abort? */
+        ckt->CKTtime += ckt->CKTdelta;
+#ifdef CLUSTER
+        CLUinput(ckt);
+#endif
+        ckt->CKTdeltaOld[0]=ckt->CKTdelta;
+        NIcomCof(ckt);
+#ifdef PREDICTOR
+        error = NIpred(ckt);
+#endif /* PREDICTOR */
+        save_mode = ckt->CKTmode;
+        save_order = ckt->CKTorder;
+#ifdef XSPICE
+/* gtri - begin - wbk - Add Breakpoint stuff */
+
+        /* Initialize temporary breakpoint to infinity */
+        g_mif_info.breakpoint.current = 1.0e30;
+
+/* gtri - end - wbk - Add Breakpoint stuff */
+
+
+/* gtri - begin - wbk - add convergence problem reporting flags */
+        /* delta is forced to equal delmin on last attempt near line 650 */
+        if(ckt->CKTdelta <= ckt->CKTdelmin)
+            ckt->enh->conv_debug.last_NIiter_call = MIF_TRUE;
+        else
+            ckt->enh->conv_debug.last_NIiter_call = MIF_FALSE;
+/* gtri - begin - wbk - add convergence problem reporting flags */
+
+
+/* gtri - begin - wbk - Call all hybrids */
+
+/* gtri - begin - wbk - Set evt_step */
+
+        if(ckt->evt->counts.num_insts > 0) {
+            g_mif_info.circuit.evt_step = ckt->CKTtime;
+        }
+/* gtri - end - wbk - Set evt_step */
+#endif
+
+        converged = NIiter(ckt,ckt->CKTtranMaxIter);
+
+#ifdef XSPICE
+        if(ckt->evt->counts.num_insts > 0) {
+            g_mif_info.circuit.evt_step = ckt->CKTtime;
+            EVTcall_hybrids(ckt);
+        }
+/* gtri - end - wbk - Call all hybrids */
+
+#endif
+        ckt->CKTstat->STATtimePts ++;
+        ckt->CKTmode = (ckt->CKTmode&MODEUIC)|MODETRAN | MODEINITPRED;
+        if(firsttime) {
+            for(i=0;i<ckt->CKTnumStates;i++) {
+                ckt->CKTstate2[i] = ckt->CKTstate1[i];
+                ckt->CKTstate3[i] = ckt->CKTstate1[i];
+            }
+        }
+        /* txl, cpl addition */
+        if (converged == 1111) {
+                return(converged);
+        }
+
+        if(converged != 0) {
+#ifndef CLUSTER
+#ifndef SHARED_MODULE
+            ckt->CKTtime = ckt->CKTtime -ckt->CKTdelta;
+            ckt->CKTstat->STATrejected ++;
+#else
+            redostep = 1;
+#endif
+#endif
+            ckt->CKTdelta = ckt->CKTdelta/8;
+#ifdef STEPDEBUG
+            (void)printf("delta cut to %g for non-convergance\n",ckt->CKTdelta);
+            fflush(stdout);
+#endif
+            if(firsttime) {
+                ckt->CKTmode = (ckt->CKTmode&MODEUIC) | MODETRAN | MODEINITTRAN;
+            }
+            ckt->CKTorder = 1;
+
+#ifdef XSPICE
+/* gtri - begin - wbk - Add Breakpoint stuff */
+
+        /* Force backup if temporary breakpoint is < current time */
+        } else if(g_mif_info.breakpoint.current < ckt->CKTtime) {
+            ckt->CKTsaveDelta = ckt->CKTdelta;
+            ckt->CKTtime -= ckt->CKTdelta;
+            ckt->CKTdelta = g_mif_info.breakpoint.current - ckt->CKTtime;
+            g_mif_info.breakpoint.last = ckt->CKTtime + ckt->CKTdelta;
+
+            if(firsttime) {
+                ckt->CKTmode = (ckt->CKTmode&MODEUIC)|MODETRAN | MODEINITTRAN;
+            }
+            ckt->CKTorder = 1;
+
+/* gtri - end - wbk - Add Breakpoint stuff */
+#endif
+
+        } else {
+
+/** Nel RHSold, ogni device deve accedere ai propri valori per vedere se esso stesso è acceso o spento.
+    Nel caso del BSIM4, vale la regola Vgs > Vth, dove Vgs = ckt->CKTrhsOld [here->BSIM4...] - ckt->CKTrhsOld [here->BSIM4...] e Vth = here->BSIM4vth .
+    In caso il transistor sia acceso, si alza un flag, privato del device, che indica che è acceso. Se è spento, lo stesso flag sarà basso.
+    Il tempo corrente CKTtime deve essere memorizzato insieme, in modo tale da poter poi calcolare il delta di tempo necessario al modello.
+    QUI, deve essere controllato che all'istante precedente il device sia acceso (o spento). Se si manifesta un cambio, allora la fase di stress (o di recovery)
+      è finita e bisogna calcolare il delta_vth attraverso il modello.
+*/
+            CKTreliability (ckt, 0) ;
+
+            if (firsttime) {
+#ifdef WANT_SENSE2
+                if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
+                    save1 = ckt->CKTmode;
+                    save2 = ckt->CKTorder;
+                    ckt->CKTmode = save_mode;
+                    ckt->CKTorder = save_order;
+                    error = CKTsenDCtran (ckt);
+                    if (error)
+                        return(error);
+
+                    ckt->CKTmode = save1;
+                    ckt->CKTorder = save2;
+                }
+#endif
+                firsttime = 0;
+#if !defined CLUSTER && !defined SHARED_MODULE
+                goto nextTime;  /* no check on
+                                 * first time point
+                                 */
+#else
+                redostep = 0;
+                goto chkStep;
+#endif
+            }
+            newdelta = ckt->CKTdelta;
+            error = CKTtrunc(ckt,&newdelta);
+            if(error) {
+                UPDATE_STATS(DOING_TRAN);
+                return(error);
+            }
+            if(newdelta > .9 * ckt->CKTdelta) {
+                if((ckt->CKTorder == 1) && (ckt->CKTmaxOrder > 1)) { /* don't rise the order for backward Euler */
+                    newdelta = ckt->CKTdelta;
+                    ckt->CKTorder = 2;
+                    error = CKTtrunc(ckt,&newdelta);
+                    if(error) {
+                        UPDATE_STATS(DOING_TRAN);
+                        return(error);
+                    }
+                    if(newdelta <= 1.05 * ckt->CKTdelta) {
+                        ckt->CKTorder = 1;
+                    }
+                }
+                /* time point OK  - 630 */
+                ckt->CKTdelta = newdelta;
+#ifdef NDEV
+                /* show a time process indicator, by Gong Ding, gdiso@ustc.edu */
+                if(ckt->CKTtime/ckt->CKTfinalTime*100<10.0)
+                    printf("%%%3.2lf\b\b\b\b\b",ckt->CKTtime/ckt->CKTfinalTime*100);
+                else  if(ckt->CKTtime/ckt->CKTfinalTime*100<100.0)
+                    printf("%%%4.2lf\b\b\b\b\b\b",ckt->CKTtime/ckt->CKTfinalTime*100);
+                else
+                    printf("%%%5.2lf\b\b\b\b\b\b\b",ckt->CKTtime/ckt->CKTfinalTime*100);
+                fflush(stdout);
+#endif
+
+#ifdef STEPDEBUG
+                (void)printf(
+                  "delta set to truncation error result: %g. Point accepted at CKTtime: %g\n",
+                  ckt->CKTdelta,ckt->CKTtime);
+                fflush(stdout);
+#endif
+
+#ifdef WANT_SENSE2
+                if(ckt->CKTsenInfo && (ckt->CKTsenInfo->SENmode & TRANSEN)){
+                    save1 = ckt->CKTmode;
+                    save2 = ckt->CKTorder;
+                    ckt->CKTmode = save_mode;
+                    ckt->CKTorder = save_order;
+                    error = CKTsenDCtran(ckt);
+                    if (error)
+                        return (error);
+
+                    ckt->CKTmode = save1;
+                    ckt->CKTorder = save2;
+                }
+#endif
+
+#if !defined CLUSTER && !defined SHARED_MODULE
+                /* go to 650 - trapezoidal */
+                goto nextTime;
+#else
+                redostep = 0;
+                goto chkStep;
+#endif
+            } else {
+#ifndef CLUSTER
+#ifndef SHARED_MODULE
+                ckt->CKTtime = ckt->CKTtime -ckt->CKTdelta;
+                ckt->CKTstat->STATrejected ++;
+#else
+                redostep = 1;
+#endif
+#endif
+                ckt->CKTdelta = newdelta;
+#ifdef STEPDEBUG
+                (void)printf(
+                    "delta set to truncation error result:point rejected\n");
+#endif
+            }
+        }
+
+        if (ckt->CKTdelta <= ckt->CKTdelmin) {
+            if (olddelta > ckt->CKTdelmin) {
+                ckt->CKTdelta = ckt->CKTdelmin;
+#ifdef STEPDEBUG
+                (void)printf("delta at delmin\n");
+#endif
+            } else {
+                UPDATE_STATS(DOING_TRAN);
+                errMsg = CKTtrouble(ckt, "Timestep too small");
+                return(E_TIMESTEP);
+            }
+        }
+#ifdef XSPICE
+/* gtri - begin - wbk - Do event backup */
+
+        if(ckt->evt->counts.num_insts > 0)
+            EVTbackup(ckt, ckt->CKTtime + ckt->CKTdelta);
+
+/* gtri - end - wbk - Do event backup */
+#endif
+#ifdef CLUSTER
+        chkStep:
+        if(CLUsync(ckt->CKTtime,&ckt->CKTdelta,redostep)){
+            goto nextTime;
+        } else {
+            ckt->CKTtime -= olddelta;
+            ckt->CKTstat->STATrejected ++;
+        }
+#endif
+
+#ifdef SHARED_MODULE
+        /* redostep == 0:
+           Either directly go to next time step, or modify ckt->CKTdelta depending on
+           synchronization requirements. sharedsync() returns 0.
+           redostep == 1:
+           No convergence, or too large truncation error.
+           Redo the last time step by subtracting olddelta, and modify ckt->CKTdelta
+           depending on synchronization requirements. sharedsync() returns 1.
+           User-supplied redo request:
+           sharedsync() may return 1 if the user has decided to do so in the callback
+           function.
+        */
+chkStep:
+        if(sharedsync(&ckt->CKTtime, &ckt->CKTdelta, olddelta, ckt->CKTfinalTime,
+                 ckt->CKTdelmin, redostep, &ckt->CKTstat->STATrejected, 1) == 0)
+            goto nextTime;
+#endif
+
+    }
+    /* NOTREACHED */
+}
diff --git a/src/spicelib/analysis/relanaskq.c b/src/spicelib/analysis/relanaskq.c
new file mode 100644
index 000000000..bb096a608
--- /dev/null
+++ b/src/spicelib/analysis/relanaskq.c
@@ -0,0 +1,44 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/ifsim.h"
+#include "ngspice/iferrmsg.h"
+#include "ngspice/relandefs.h"
+#include "ngspice/cktdefs.h"
+
+int 
+RELANaskQuest (CKTcircuit *ckt, JOB *anal, int which, IFvalue *value)
+{
+    RELANan *job = (RELANan *) anal ;
+
+    NG_IGNORE (ckt) ;
+
+    switch (which)
+    {
+        case RELAN_TSTART:
+            value->rValue = job->RELANinitTime ;
+            break ;
+        case RELAN_TSTOP:
+            value->rValue = job->RELANfinalTime ;
+            break ;
+        case RELAN_TSTEP:
+            value->rValue = job->RELANstep ;
+            break ;
+        case RELAN_TMAX:
+            value->rValue = job->RELANmaxStep ;
+            break ;
+        case RELAN_UIC:
+            if (job->RELANmode & MODEUIC) {
+                value->iValue = 1 ;
+            } else {
+                value->iValue = 0 ;
+            }
+            break ;
+        default:
+            return (E_BADPARM) ;
+    }
+
+    return (OK) ;
+}
diff --git a/src/spicelib/analysis/relaninit.c b/src/spicelib/analysis/relaninit.c
new file mode 100644
index 000000000..9b4b4c759
--- /dev/null
+++ b/src/spicelib/analysis/relaninit.c
@@ -0,0 +1,34 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/cktdefs.h"
+#include "ngspice/relandefs.h"
+#include "ngspice/iferrmsg.h"
+
+int
+RELANinit (CKTcircuit *ckt, JOB *anal)
+{
+    RELANan *job = (RELANan *) anal ;
+
+    ckt->CKTinitTime  = job->RELANinitTime ;
+    ckt->CKTfinalTime = job->RELANfinalTime ;
+    ckt->CKTstep      = job->RELANstep ;
+    ckt->CKTmaxStep   = job->RELANmaxStep ;
+
+    if (ckt->CKTmaxStep == 0)
+    {
+        if (ckt->CKTstep < (ckt->CKTfinalTime - ckt->CKTinitTime) / 50.0)
+        {
+            ckt->CKTmaxStep = ckt->CKTstep ;
+        } else {
+            ckt->CKTmaxStep = (ckt->CKTfinalTime - ckt->CKTinitTime) / 50.0 ;
+        }
+    }
+
+    ckt->CKTdelmin = 1e-11 * ckt->CKTmaxStep ;
+    ckt->CKTmode = job->RELANmode ;
+
+    return (OK) ;
+}
diff --git a/src/spicelib/analysis/relansetp.c b/src/spicelib/analysis/relansetp.c
new file mode 100644
index 000000000..73c67623c
--- /dev/null
+++ b/src/spicelib/analysis/relansetp.c
@@ -0,0 +1,89 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/ifsim.h"
+#include "ngspice/iferrmsg.h"
+#include "ngspice/relandefs.h"
+#include "ngspice/cktdefs.h"
+
+#include "analysis.h"
+
+int 
+RELANsetParm (CKTcircuit *ckt, JOB *anal, int which, IFvalue *value)
+{
+    RELANan *job = (RELANan *) anal ;
+
+    NG_IGNORE (ckt) ;
+
+    switch (which)
+    {
+        case RELAN_TSTART:
+            if (value->rValue >= job->RELANfinalTime)
+            {
+                errMsg = copy ("TSTART is invalid, must be less than TSTOP.") ;
+                job->RELANinitTime = 0.0 ;
+                return (E_PARMVAL) ;
+            }
+            job->RELANinitTime = value->rValue ;
+            break ;
+        case RELAN_TSTOP:
+            if (value->rValue <= 0.0)
+            {
+                errMsg = copy ("TST0P is invalid, must be greater than zero.") ;
+                job->RELANfinalTime = 1.0 ;
+                return (E_PARMVAL) ;
+            }
+            job->RELANfinalTime = value->rValue ;
+            break ;
+        case RELAN_TSTEP:
+            if (value->rValue <= 0.0)
+            {
+                errMsg = copy ("TSTEP is invalid, must be greater than zero.") ;
+                job->RELANstep = 1.0 ;
+                return (E_PARMVAL) ;
+	    }
+            job->RELANstep = value->rValue ;
+            break ;
+        case RELAN_TMAX:
+            job->RELANmaxStep = value->rValue ;
+            break ;
+        case RELAN_UIC:
+            if (value->iValue)
+            {
+                job->RELANmode |= MODEUIC ;
+            }
+            break ;
+        default:
+            return (E_BADPARM) ;
+    }
+
+    return (OK) ;
+}
+
+
+static IFparm RELANparms [] = {
+    { "relan_aging_start",     RELAN_TSTART,    IF_SET|IF_REAL, "starting time" },
+    { "relan_aging_stop",      RELAN_TSTOP,     IF_SET|IF_REAL, "ending time" },
+    { "relan_aging_step",      RELAN_TSTEP,     IF_SET|IF_REAL, "time step" },
+    { "relan_aging_max",       RELAN_TMAX,      IF_SET|IF_REAL, "maximum time step" },
+    { "uic",                   RELAN_UIC,       IF_SET|IF_FLAG, "use initial conditions" },
+} ;
+
+SPICEanalysis RELANinfo = {
+    { 
+        "RELAN",
+        "Reliability analysis",
+
+        NUMELEMS(RELANparms),
+        RELANparms
+    },
+    sizeof(RELANan),
+    TIMEDOMAIN,
+    1,
+    RELANsetParm,
+    RELANaskQuest,
+    RELANinit,
+    RELANanalysis
+} ;
diff --git a/src/spicelib/devices/bsim4/Makefile.am b/src/spicelib/devices/bsim4/Makefile.am
index 047ffdfc6..483dccb86 100644
--- a/src/spicelib/devices/bsim4/Makefile.am
+++ b/src/spicelib/devices/bsim4/Makefile.am
@@ -27,7 +27,8 @@ libbsim4_la_SOURCES =	\
 	bsim4ext.h	\
 	bsim4init.c	\
 	bsim4init.h	\
-	bsim4itf.h
+	bsim4itf.h	\
+	b4reliability.c
 
 
 
diff --git a/src/spicelib/devices/bsim4/b4reliability.c b/src/spicelib/devices/bsim4/b4reliability.c
new file mode 100644
index 000000000..65d2386cb
--- /dev/null
+++ b/src/spicelib/devices/bsim4/b4reliability.c
@@ -0,0 +1,241 @@
+/**********
+Author: 2015 Francesco Lannutti
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/devdefs.h"
+#include "bsim4def.h"
+#include "ngspice/sperror.h"
+
+static int
+calculate_aging
+(
+    BSIM4instance *here,
+    double t_aging,
+    unsigned int stress_or_recovery
+)
+{
+    double K_b, T, h_cut, q, Nts, T_hk, Nt, eps_0, eps_hk, eps_SiO2, m_star, W, tau_0, beta, tau_e, beta1 ;
+
+    double A ;
+
+    K_b = 8.6e-5 ;
+    T = 300 ;
+    h_cut = 1.05e-34 ;
+    q = 1.6e-19 ;
+    Nts = 2e13 ;
+    T_hk = 2 ;
+    Nt = pow ((sqrt (Nts)), 3) * 1e-21 ;
+    eps_0 = 8.85e-21 ;
+    eps_hk = 25 ;
+    eps_SiO2 = 3.9 ;
+    m_star = 0.1 * 9.11e-31 ;
+    W = 1.5 * 1.6e-19 ;
+    tau_0 = 1e-11 ;
+    beta = 0.373 ;
+    tau_e = 0.85e-9 ;
+    beta1 = 0.112 ;
+
+    A = (q / (4 * eps_0 * eps_hk)) * pow ((h_cut / (2 * sqrt (2 * m_star * W)) * 1e9), 2) ;
+
+    if (stress_or_recovery)
+    {
+        if (h_cut / (2 * sqrt (2 * m_star * W)) * log (1 + pow (((t_aging + here->BSIM4reliability->t_star) / tau_0), beta)) * 1e9 <= 2)
+        {
+            here->BSIM4reliability->deltaVth = Nt * A * pow (log (1 + pow (((t_aging + here->BSIM4reliability->t_star) / tau_0), beta)), 2) ;
+        } else {
+            here->BSIM4reliability->deltaVth = pow ((q / (4 * eps_0 * eps_hk)) * Nt * T_hk, 2) ;
+        }
+    } else {
+        here->BSIM4reliability->deltaVth = here->BSIM4reliability->deltaVth * log (1 + (1.718 / (1 + pow ((t_aging / tau_e), beta1)))) ;
+    }
+
+    if (!stress_or_recovery)
+    {
+        here->BSIM4reliability->t_star = pow ((exp (sqrt (here->BSIM4reliability->deltaVth / (Nt * A))) - 1), (1 / beta)) * tau_0 ;
+    }
+
+    return 0 ;
+}
+
+int
+BSIM4reliability (GENmodel *inModel, CKTcircuit *ckt, unsigned int mode)
+{
+    BSIM4model *model = (BSIM4model *)inModel ;
+    BSIM4instance *here ;
+    double delta, vds, vgs, von ;
+    int NowIsON ;
+
+    /*  loop through all the BSIM4 device models */
+    for ( ; model != NULL ; model = model->BSIM4nextModel)
+    {
+        /* loop through all the instances of the model */
+        for (here = model->BSIM4instances ; here != NULL ; here=here->BSIM4nextInstance)
+        {
+            if (model->BSIM4type == -1)
+            {
+                // Determine if the transistor is ON or OFF
+                vds = ckt->CKTstate0 [here->BSIM4vds] ;
+                vgs = ckt->CKTstate0 [here->BSIM4vgs] ;
+                von = model->BSIM4type * here->BSIM4von ;
+                if (vds >= 0)
+                {
+//                    printf ("VDS >= 0\tBSIM4type: %d\tBSIM4instance: %s\tVgs: %-.9g\tVon: %-.9g\t", model->BSIM4type, here->BSIM4name, vgs, von) ;
+                    if (vgs > von)
+                    {
+                        if (here->BSIM4rgateMod == 3)
+                        {
+                            double vges, vgms ;
+                            vges = ckt->CKTstate0 [here->BSIM4vges] ;
+                            vgms = ckt->CKTstate0 [here->BSIM4vgms] ;
+                            if ((vges > von) && (vgms > von))
+                            {
+//                                printf ("Acceso!\n") ;
+                                NowIsON = 1 ;
+                            } else {
+//                                printf ("Spento!\n") ;
+                                NowIsON = 0 ;
+                            }
+                        } else if ((here->BSIM4rgateMod == 1) || (here->BSIM4rgateMod == 2)) {
+                            double vges ;
+                            vges = ckt->CKTstate0 [here->BSIM4vges] ;
+                            if (vges > von)
+                            {
+//                                printf ("Acceso!\n") ;
+                                NowIsON = 1 ;
+                            } else {
+//                                printf ("Spento!\n") ;
+                                NowIsON = 0 ;
+                            }
+                        } else {
+//                            printf ("Acceso!\n") ;
+                            NowIsON = 1 ;
+                        }
+                    } else {
+//                        printf ("Spento!\n") ;
+                        NowIsON = 0 ;
+                    }
+                } else {
+                    double vgd ;
+                    vgd = vgs - vds ;
+//                    printf ("VDS <  0\tBSIM4type: %d\tBSIM4instance: %s\tVgd: %-.9g\tVon: %-.9g\t", model->BSIM4type, here->BSIM4name, vgd, von) ;
+                    if (vgd > von)
+                    {
+                        if (here->BSIM4rgateMod == 3)
+                        {
+                            double vges, vged, vgms, vgmd ;
+                            vges = ckt->CKTstate0 [here->BSIM4vges] ;
+                            vged = vges - vds ;
+                            vgms = ckt->CKTstate0 [here->BSIM4vgms] ;
+                            vgmd = vgms - vds ;
+                            if ((vged > von) && (vgmd > von))
+                            {
+//                                printf ("Acceso!\n") ;
+                                NowIsON = 1 ;
+                            } else {
+//                                printf ("Spento!\n") ;
+                                NowIsON = 0 ;
+                            }
+                        } else if ((here->BSIM4rgateMod == 1) || (here->BSIM4rgateMod == 2)) {
+                            double vges, vged ;
+                            vges = ckt->CKTstate0 [here->BSIM4vges] ;
+                            vged = vges - vds ;
+                            if (vged > von)
+                            {
+//                                printf ("Acceso!\n") ;
+                                NowIsON = 1 ;
+                            } else {
+//                                printf ("Spento!\n") ;
+                                NowIsON = 0 ;
+                            }
+                        } else {
+//                            printf ("Acceso!\n") ;
+                            NowIsON = 1 ;
+                        }
+                    } else {
+//                        printf ("Spento!\n") ;
+                        NowIsON = 0 ;
+                    }
+                }
+
+                // If it's the first time, initialize 'here->BSIM4reliability->IsON'
+                if (here->BSIM4reliability->IsON == -1)
+                {
+                    here->BSIM4reliability->IsON = NowIsON ;
+                }
+
+                if (mode == 0)
+                {
+                    if (NowIsON)
+                    {
+                        if (here->BSIM4reliability->IsON == 1)
+                        {
+                            // Until now, the device was ON - Do NOTHING
+                            delta = -1 ;
+                        } else if (here->BSIM4reliability->IsON == 0) {
+                            // Until now, the device was OFF - Calculate recovery
+                            delta = ckt->CKTtime - here->BSIM4reliability->time ;
+
+                            // Calculate Aging - Giogio Liatis' Model
+                            calculate_aging (here, delta, 0) ;
+
+                            // Update time and flag - Stress begins
+                            here->BSIM4reliability->time = ckt->CKTtime ;
+                            here->BSIM4reliability->IsON = 1 ;
+                        } else {
+                            fprintf (stderr, "Reliability Analysis Error\n") ;
+                        }
+                    } else {
+                        if (here->BSIM4reliability->IsON == 1)
+                        {
+                            // Until now, the device was ON - Calculate stress
+                            delta = ckt->CKTtime - here->BSIM4reliability->time ;
+
+                            // Calculate Aging - Giorgio Liatis' Model
+                            calculate_aging (here, delta, 1) ;
+
+                            // Update time and flag - Recovery begins
+                            here->BSIM4reliability->time = ckt->CKTtime ;
+                            here->BSIM4reliability->IsON = 0 ;
+                        } else if (here->BSIM4reliability->IsON == 0) {
+                            // Until now, the device was OFF - Do NOTHING
+                            delta = -1 ;
+                        } else {
+                            fprintf (stderr, "Reliability Analysis Error\n") ;
+                        }
+                    }
+                } else if (mode == 1) {
+                    // In this mode, it doesn't matter if NOW the device is in stress or in recovery, since it's the last timestep
+                    if (here->BSIM4reliability->IsON == 1)
+                    {
+                        // Calculate stress
+                        delta = ckt->CKTtime - here->BSIM4reliability->time ;
+                        calculate_aging (here, delta, 1) ;
+
+                        // Update time and flag - Maybe Optional
+                        here->BSIM4reliability->time = ckt->CKTtime ;
+                        here->BSIM4reliability->IsON = 1 ;
+                    } else if (here->BSIM4reliability->IsON == 0) {
+                        // Calculate recovery
+                        delta = ckt->CKTtime - here->BSIM4reliability->time ;
+                        calculate_aging (here, delta, 0) ;
+
+                        // Update time and flag - Maybe Optional
+                        here->BSIM4reliability->time = ckt->CKTtime ;
+                        here->BSIM4reliability->IsON = 0 ;
+                    } else {
+                        fprintf (stderr, "Reliability Analysis Error\n") ;
+                    }
+                    printf ("\tTime: %-.9gs\t\t", ckt->CKTtime) ;
+                    printf ("DeltaVth: %-.9gmV\t\t", here->BSIM4reliability->deltaVth * 1000) ;
+                    printf ("Device Name: %s\t\t", here->BSIM4name) ;
+                    printf ("Device Type: %s\n\n", model->BSIM4modName) ;
+                } else {
+                    fprintf (stderr, "Reliability Analysis Error\n") ;
+                }
+            }
+        }
+    }
+
+    return (OK) ;
+}
diff --git a/src/spicelib/devices/bsim4/b4set.c b/src/spicelib/devices/bsim4/b4set.c
index fa445937d..1f8349531 100644
--- a/src/spicelib/devices/bsim4/b4set.c
+++ b/src/spicelib/devices/bsim4/b4set.c
@@ -2653,6 +2653,15 @@ do { if((here->ptr = SMPmakeElt(matrix,here->first,here->second))==(double *)NUL
                 TSTALLOC(BSIM4SgpPtr, BSIM4sNode, BSIM4gNodePrime);
                 TSTALLOC(BSIM4SbpPtr, BSIM4sNode, BSIM4bNodePrime);
             }
+
+            if (model->BSIM4type == -1)
+            {
+                here->BSIM4reliability = TMALLOC (BSIM4relStruct, 1) ;
+                here->BSIM4reliability->time = 0 ;
+                here->BSIM4reliability->deltaVth = 0 ;
+                here->BSIM4reliability->t_star = 0 ;
+                here->BSIM4reliability->IsON = -1 ;
+            }
         }
     }
 
diff --git a/src/spicelib/devices/bsim4/b4temp.c b/src/spicelib/devices/bsim4/b4temp.c
index 9b71b127d..a446e2032 100644
--- a/src/spicelib/devices/bsim4/b4temp.c
+++ b/src/spicelib/devices/bsim4/b4temp.c
@@ -1715,6 +1715,13 @@ int Size_Not_Found, i;
               }
 
             /* adding delvto  */
+            if (model->BSIM4type == -1)
+            {
+                if (here->BSIM4reliability->IsON != -1)
+                {
+                    here->BSIM4delvto = here->BSIM4reliability->deltaVth ;
+                }
+            }
             here->BSIM4vth0 += here->BSIM4delvto;
             here->BSIM4vfb = pParam->BSIM4vfb + model->BSIM4type * here->BSIM4delvto;
 
diff --git a/src/spicelib/devices/bsim4/bsim4def.h b/src/spicelib/devices/bsim4/bsim4def.h
index 0b4efec63..59e2f31d9 100644
--- a/src/spicelib/devices/bsim4/bsim4def.h
+++ b/src/spicelib/devices/bsim4/bsim4def.h
@@ -29,6 +29,13 @@ Modified by Pankaj Kumar Thakur, 07/23/2012
 #include "ngspice/complex.h"
 #include "ngspice/noisedef.h"
 
+typedef struct sBSIM4relStruct {
+    double time ;
+    double deltaVth ;
+    double t_star ;
+    int IsON ;
+} BSIM4relStruct ;
+
 typedef struct sBSIM4instance
 {
     struct sBSIM4model *BSIM4modPtr;
@@ -574,6 +581,8 @@ typedef struct sBSIM4instance
         double **BSIM4nVar;
 #endif /* NONOISE */
 
+    BSIM4relStruct *BSIM4reliability ;
+
 } BSIM4instance ;
 
 struct bsim4SizeDependParam
diff --git a/src/spicelib/devices/bsim4/bsim4ext.h b/src/spicelib/devices/bsim4/bsim4ext.h
index 9618c6328..693778691 100644
--- a/src/spicelib/devices/bsim4/bsim4ext.h
+++ b/src/spicelib/devices/bsim4/bsim4ext.h
@@ -31,3 +31,5 @@ extern int BSIM4trunc(GENmodel*,CKTcircuit*,double*);
 extern int BSIM4noise(int,int,GENmodel*,CKTcircuit*,Ndata*,double*);
 extern int BSIM4unsetup(GENmodel*,CKTcircuit*);
 extern int BSIM4soaCheck(CKTcircuit *, GENmodel *);
+
+extern int BSIM4reliability (GENmodel *, CKTcircuit *, unsigned int) ;
diff --git a/src/spicelib/devices/bsim4/bsim4init.c b/src/spicelib/devices/bsim4/bsim4init.c
index e9d13577d..b355a6641 100644
--- a/src/spicelib/devices/bsim4/bsim4init.c
+++ b/src/spicelib/devices/bsim4/bsim4init.c
@@ -73,7 +73,10 @@ SPICEdev BSIM4info = {
     NULL,          /* DEVacct        */
 #endif
     &BSIM4iSize,   /* DEVinstSize    */
-    &BSIM4mSize    /* DEVmodSize     */
+    &BSIM4mSize,    /* DEVmodSize     */
+
+ BSIM4reliability
+
 };
 
 
diff --git a/src/spicelib/devices/mos1/mos1reliability.c b/src/spicelib/devices/mos1/mos1reliability.c
index 2068f34b9..103d6b055 100644
--- a/src/spicelib/devices/mos1/mos1reliability.c
+++ b/src/spicelib/devices/mos1/mos1reliability.c
@@ -79,25 +79,25 @@ MOS1reliability (GENmodel *inModel, CKTcircuit *ckt, unsigned int mode)
                 von = model->MOS1type * here->MOS1von ;
                 if (vds >= 0)
                 {
-                    printf ("VDS >= 0\tMOS1type: %d\tMOS1instance: %s\tVgs: %-.9g\tVon: %-.9g\t", model->MOS1type, here->MOS1name, vgs, von) ;
+//                    printf ("VDS >= 0\tMOS1type: %d\tMOS1instance: %s\tVgs: %-.9g\tVon: %-.9g\t", model->MOS1type, here->MOS1name, vgs, von) ;
                     if (vgs > von)
                     {
-                        printf ("Acceso!\n") ;
+//                        printf ("Acceso!\n") ;
                         NowIsON = 1 ;
                     } else {
-                        printf ("Spento!\n") ;
+//                        printf ("Spento!\n") ;
                         NowIsON = 0 ;
                     }
                 } else {
                     double vgd ;
                     vgd = vgs - vds ;
-                    printf ("VDS <  0\tMOS1type: %d\tMOS1instance: %s\tVgd: %-.9g\tVon: %-.9g\t", model->MOS1type, here->MOS1name, vgd, von) ;
+//                    printf ("VDS <  0\tMOS1type: %d\tMOS1instance: %s\tVgd: %-.9g\tVon: %-.9g\t", model->MOS1type, here->MOS1name, vgd, von) ;
                     if (vgd > von)
                     {
-                        printf ("Acceso!\n") ;
+//                        printf ("Acceso!\n") ;
                         NowIsON = 1 ;
                     } else {
-                        printf ("Spento!\n") ;
+//                        printf ("Spento!\n") ;
                         NowIsON = 0 ;
                     }
                 }
@@ -170,13 +170,13 @@ MOS1reliability (GENmodel *inModel, CKTcircuit *ckt, unsigned int mode)
                     } else {
                         fprintf (stderr, "Reliability Analysis Error\n") ;
                     }
+                    printf ("\tTime: %-.9gs\t\t", ckt->CKTtime) ;
+                    printf ("DeltaVth: %-.9gmV\t\t", here->MOS1reliability->deltaVth * 1000) ;
+                    printf ("Device Name: %s\t\t", here->MOS1name) ;
+                    printf ("Device Type: %s\n\n", model->MOS1modName) ;
                 } else {
                     fprintf (stderr, "Reliability Analysis Error\n") ;
                 }
-                printf ("Time: %-.9gs\t", here->MOS1reliability->time) ;
-                printf ("DeltaVth: %-.9gmV\t", here->MOS1reliability->deltaVth * 1000) ;
-                printf ("IsON: %u\t", here->MOS1reliability->IsON) ;
-                printf ("t_star: %-.9gs\n\n\n", here->MOS1reliability->t_star) ;
             }
         }
     }
diff --git a/src/spicelib/parser/inp2dot.c b/src/spicelib/parser/inp2dot.c
index fdc64e2a3..0f0872e71 100644
--- a/src/spicelib/parser/inp2dot.c
+++ b/src/spicelib/parser/inp2dot.c
@@ -612,6 +612,66 @@ dot_pss(char *line, void *ckt, INPtables *tab, card *current,
 /* SP */
 #endif
 
+//#ifdef RELAN
+static int
+dot_relan (char *line, void *ckt, INPtables *tab, card *current, void *task, void *gnode, JOB *foo)
+{
+    int error ;      /* error code temporary */
+    IFvalue ptemp ;  /* a value structure to package resistance into */
+    IFvalue *parm ;  /* a pointer to a value struct for function returns */
+    int which ;      /* which analysis we are performing */
+    char *word ;     /* something to stick a word of input into */
+    double dtemp ;   /* random double precision temporary */
+
+    NG_IGNORE (gnode) ;
+
+    /* .relan AgingStep AgingStop <AgingStart> <tmax> <UIC> */
+    which = ft_find_analysis ("RELAN") ;
+    if (which == -1)
+    {
+        LITERR ("Reliability Analysis unsupported.\n") ;
+        return (0) ;
+    }
+    IFC (newAnalysis, (ckt, which, "Reliability Analysis", &foo, task)) ;
+
+    parm = INPgetValue (ckt, &line, IF_REAL, tab) ;
+    GCA (INPapName, (ckt, which, foo, "relan_aging_step", parm)) ;
+
+    parm = INPgetValue (ckt, &line, IF_REAL, tab) ;
+    GCA (INPapName, (ckt, which, foo, "relan_aging_stop", parm)) ;
+
+    if (*line)
+    {
+        dtemp = INPevaluate (&line, &error, 1) ; /* AgingStart */
+        if (error == 0)
+        {
+            ptemp.rValue = dtemp ;
+            GCA (INPapName, (ckt, which, foo, "relan_aging_start", &ptemp)) ;
+            dtemp = INPevaluate (&line, &error, 1) ;	/* tmax */
+            if (error == 0)
+            {
+                ptemp.rValue = dtemp ;
+                GCA (INPapName, (ckt, which, foo, "relan_aging_max", &ptemp)) ;
+            }
+        }
+    }
+
+    if (*line)
+    {
+        INPgetTok (&line, &word, 1) ; /* UIC */
+        if (strcmp (word, "uic") == 0)
+        {
+            ptemp.iValue = 1 ;
+            GCA (INPapName, (ckt, which, foo, "uic", &ptemp)) ;
+        } else {
+            fprintf (stderr,"Error: unknown parameter %s on .relan - ignored\n", word) ;
+        }
+    }
+
+    return (0) ;
+}
+//#endif
+
 static int
 dot_options(char *line, CKTcircuit *ckt, INPtables *tab, card *current,
             TSKtask *task, CKTnode *gnode, JOB *foo)
@@ -698,6 +758,13 @@ INP2dot(CKTcircuit *ckt, INPtables *tab, card *current, TSKtask *task, CKTnode *
 		goto quit;
 		/* SP */
 #endif
+
+//#ifdef RELAN
+        } else if ((strcmp (token, ".relan") == 0)) {
+		rtn = dot_relan (line, ckt, tab, current, task, gnode, foo) ;
+		goto quit ;
+//#endif
+
 	} else if ((strcmp(token, ".subckt") == 0) ||
 	           (strcmp(token, ".ends") == 0)) {
 		/* not yet implemented - warn & ignore */