diff --git a/configure.ac b/configure.ac
index 47c7d93f2..de5ba22ad 100644
--- a/configure.ac
+++ b/configure.ac
@@ -1112,6 +1112,7 @@ AC_CONFIG_FILES([Makefile
                  src/spicelib/devices/mos6/Makefile
                  src/spicelib/devices/mos9/Makefile
                  src/spicelib/devices/ndev/Makefile
+                 src/spicelib/devices/relmodel/Makefile
                  src/spicelib/devices/res/Makefile
                  src/spicelib/devices/soi3/Makefile
                  src/spicelib/devices/sw/Makefile
diff --git a/src/Makefile.am b/src/Makefile.am
index eec46cd59..17c884461 100644
--- a/src/Makefile.am
+++ b/src/Makefile.am
@@ -138,6 +138,11 @@ ngspice_LDADD += \
 	spicelib/devices/numos/libnumos.la
 endif
 
+if RELAN_WANTED
+ngspice_LDADD += \
+	spicelib/devices/relmodel/librelmodel.la
+endif
+
 ngspice_LDADD += \
 	spicelib/analysis/libckt.la \
 	spicelib/devices/libdev.la
@@ -385,6 +390,11 @@ libspice_la_LIBADD += \
 	spicelib/devices/numos/libnumos.la
 endif
 
+if RELAN_WANTED
+libspice_la_LIBADD += \
+	spicelib/devices/relmodel/librelmodel.la
+endif
+
 libspice_la_LIBADD += \
 	spicelib/analysis/libckt.la \
 	spicelib/devices/libdev.la
@@ -500,6 +510,11 @@ libngspice_la_LIBADD += \
 	spicelib/devices/numos/libnumos.la
 endif
 
+if RELAN_WANTED
+libngspice_la_LIBADD += \
+	spicelib/devices/relmodel/librelmodel.la
+endif
+
 libngspice_la_LIBADD += \
 	spicelib/analysis/libckt.la \
 	spicelib/devices/libdev.la
diff --git a/src/frontend/inpcom.c b/src/frontend/inpcom.c
index 69bea46dd..f6384a018 100644
--- a/src/frontend/inpcom.c
+++ b/src/frontend/inpcom.c
@@ -1877,6 +1877,11 @@ comment_out_unused_subckt_models(struct line *start_card)
             if (!cieq(model_type, "c") &&
                 !cieq(model_type, "l") &&
                 !cieq(model_type, "r") &&
+
+#ifdef RELAN
+                !cieq(model_type, "relmodel") &&
+#endif
+
                 !nlist_model_find(used_models, model_name))
             {
                 *line = '*';
diff --git a/src/include/ngspice/gendefs.h b/src/include/ngspice/gendefs.h
index 8bc1546cc..0546d22f0 100644
--- a/src/include/ngspice/gendefs.h
+++ b/src/include/ngspice/gendefs.h
@@ -40,6 +40,11 @@ struct GENmodel {       /* model structure for a resistor */
     GENinstance *GENinstances;  /* pointer to list of instances that have this
                                  * model */
     IFuid GENmodName;           /* pointer to character string naming this model */
+
+#ifdef RELAN
+    GENmodel *GENrelmodelModel ; /* Relmodel Companion Model */
+#endif
+
 };
 
 #endif
diff --git a/src/spicelib/devices/Makefile.am b/src/spicelib/devices/Makefile.am
index 706738dee..f4c85949b 100644
--- a/src/spicelib/devices/Makefile.am
+++ b/src/spicelib/devices/Makefile.am
@@ -60,6 +60,10 @@ if NUMDEV_WANTED
 SUBDIRS += nbjt nbjt2 numd  numd2 numos
 endif
 
+if RELAN_WANTED
+SUBDIRS += relmodel
+endif
+
 DIST_SUBDIRS = \
 	@VLADEVDIR@ \
 	asrc	\
diff --git a/src/spicelib/devices/bsim4/b4reliability.c b/src/spicelib/devices/bsim4/b4reliability.c
index 65d2386cb..afe700a1f 100644
--- a/src/spicelib/devices/bsim4/b4reliability.c
+++ b/src/spicelib/devices/bsim4/b4reliability.c
@@ -1,5 +1,5 @@
 /**********
-Author: 2015 Francesco Lannutti
+Author: 2015 Francesco Lannutti - July 2015
 **********/
 
 #include "ngspice/ngspice.h"
@@ -7,64 +7,13 @@ Author: 2015 Francesco Lannutti
 #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 ;
+    int NowIsON, ret ;
 
     /*  loop through all the BSIM4 device models */
     for ( ; model != NULL ; model = model->BSIM4nextModel)
@@ -72,7 +21,7 @@ BSIM4reliability (GENmodel *inModel, CKTcircuit *ckt, unsigned int mode)
         /* loop through all the instances of the model */
         for (here = model->BSIM4instances ; here != NULL ; here=here->BSIM4nextInstance)
         {
-            if (model->BSIM4type == -1)
+            if (model->BSIM4type == PMOS)
             {
                 // Determine if the transistor is ON or OFF
                 vds = ckt->CKTstate0 [here->BSIM4vds] ;
@@ -158,46 +107,54 @@ BSIM4reliability (GENmodel *inModel, CKTcircuit *ckt, unsigned int mode)
                     }
                 }
 
-                // If it's the first time, initialize 'here->BSIM4reliability->IsON'
-                if (here->BSIM4reliability->IsON == -1)
+                // If it's the first time, initialize 'here->relStruct->IsON'
+                if (here->relStruct->IsON == -1)
                 {
-                    here->BSIM4reliability->IsON = NowIsON ;
+                    here->relStruct->IsON = NowIsON ;
                 }
 
                 if (mode == 0)
                 {
                     if (NowIsON)
                     {
-                        if (here->BSIM4reliability->IsON == 1)
+                        if (here->relStruct->IsON == 1)
                         {
                             // Until now, the device was ON - Do NOTHING
                             delta = -1 ;
-                        } else if (here->BSIM4reliability->IsON == 0) {
+                        } else if (here->relStruct->IsON == 0) {
                             // Until now, the device was OFF - Calculate recovery
-                            delta = ckt->CKTtime - here->BSIM4reliability->time ;
+                            delta = ckt->CKTtime - here->relStruct->time ;
 
                             // Calculate Aging - Giogio Liatis' Model
-                            calculate_aging (here, delta, 0) ;
+                            ret = RELMODELcalculateAging ((GENinstance *)here, here->BSIM4modPtr->BSIM4modType, delta, 0) ;
+                            if (ret == 1)
+                            {
+                                return (E_INTERN) ;
+                            }
 
                             // Update time and flag - Stress begins
-                            here->BSIM4reliability->time = ckt->CKTtime ;
-                            here->BSIM4reliability->IsON = 1 ;
+                            here->relStruct->time = ckt->CKTtime ;
+                            here->relStruct->IsON = 1 ;
                         } else {
                             fprintf (stderr, "Reliability Analysis Error\n") ;
                         }
                     } else {
-                        if (here->BSIM4reliability->IsON == 1)
+                        if (here->relStruct->IsON == 1)
                         {
                             // Until now, the device was ON - Calculate stress
-                            delta = ckt->CKTtime - here->BSIM4reliability->time ;
+                            delta = ckt->CKTtime - here->relStruct->time ;
 
                             // Calculate Aging - Giorgio Liatis' Model
-                            calculate_aging (here, delta, 1) ;
+                            ret = RELMODELcalculateAging ((GENinstance *)here, here->BSIM4modPtr->BSIM4modType, delta, 1) ;
+                            if (ret == 1)
+                            {
+                                return (E_INTERN) ;
+                            }
 
                             // Update time and flag - Recovery begins
-                            here->BSIM4reliability->time = ckt->CKTtime ;
-                            here->BSIM4reliability->IsON = 0 ;
-                        } else if (here->BSIM4reliability->IsON == 0) {
+                            here->relStruct->time = ckt->CKTtime ;
+                            here->relStruct->IsON = 0 ;
+                        } else if (here->relStruct->IsON == 0) {
                             // Until now, the device was OFF - Do NOTHING
                             delta = -1 ;
                         } else {
@@ -206,28 +163,36 @@ BSIM4reliability (GENmodel *inModel, CKTcircuit *ckt, unsigned int mode)
                     }
                 } 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)
+                    if (here->relStruct->IsON == 1)
                     {
                         // Calculate stress
-                        delta = ckt->CKTtime - here->BSIM4reliability->time ;
-                        calculate_aging (here, delta, 1) ;
+                        delta = ckt->CKTtime - here->relStruct->time ;
+                        ret = RELMODELcalculateAging ((GENinstance *)here, here->BSIM4modPtr->BSIM4modType, delta, 1) ;
+                        if (ret == 1)
+                        {
+                            return (E_INTERN) ;
+                        }
 
                         // Update time and flag - Maybe Optional
-                        here->BSIM4reliability->time = ckt->CKTtime ;
-                        here->BSIM4reliability->IsON = 1 ;
-                    } else if (here->BSIM4reliability->IsON == 0) {
+                        here->relStruct->time = ckt->CKTtime ;
+                        here->relStruct->IsON = 1 ;
+                    } else if (here->relStruct->IsON == 0) {
                         // Calculate recovery
-                        delta = ckt->CKTtime - here->BSIM4reliability->time ;
-                        calculate_aging (here, delta, 0) ;
+                        delta = ckt->CKTtime - here->relStruct->time ;
+                        ret = RELMODELcalculateAging ((GENinstance *)here, here->BSIM4modPtr->BSIM4modType, delta, 0) ;
+                        if (ret == 1)
+                        {
+                            return (E_INTERN) ;
+                        }
 
                         // Update time and flag - Maybe Optional
-                        here->BSIM4reliability->time = ckt->CKTtime ;
-                        here->BSIM4reliability->IsON = 0 ;
+                        here->relStruct->time = ckt->CKTtime ;
+                        here->relStruct->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 ("DeltaVth: %-.9gmV\t\t", here->relStruct->deltaVth * 1000) ;
                     printf ("Device Name: %s\t\t", here->BSIM4name) ;
                     printf ("Device Type: %s\n\n", model->BSIM4modName) ;
                 } else {
diff --git a/src/spicelib/devices/bsim4/b4set.c b/src/spicelib/devices/bsim4/b4set.c
index bfb9a8a73..687aad5b3 100644
--- a/src/spicelib/devices/bsim4/b4set.c
+++ b/src/spicelib/devices/bsim4/b4set.c
@@ -2655,13 +2655,14 @@ do { if((here->ptr = SMPmakeElt(matrix,here->first,here->second))==(double *)NUL
             }
 
 #ifdef RELAN
-            if (model->BSIM4type == -1)
+            here->relStruct = NULL ;
+            if (model->BSIM4type == PMOS)
             {
-                here->BSIM4reliability = TMALLOC (BSIM4relStruct, 1) ;
-                here->BSIM4reliability->time = 0 ;
-                here->BSIM4reliability->deltaVth = 0 ;
-                here->BSIM4reliability->t_star = 0 ;
-                here->BSIM4reliability->IsON = -1 ;
+                here->relStruct = TMALLOC (RELMODELrelStruct, 1) ;
+                here->relStruct->time = 0 ;
+                here->relStruct->deltaVth = 0 ;
+                here->relStruct->t_star = 0 ;
+                here->relStruct->IsON = -1 ;
             }
 #endif
 
diff --git a/src/spicelib/devices/bsim4/b4temp.c b/src/spicelib/devices/bsim4/b4temp.c
index b7186c774..81908a030 100644
--- a/src/spicelib/devices/bsim4/b4temp.c
+++ b/src/spicelib/devices/bsim4/b4temp.c
@@ -1717,11 +1717,11 @@ int Size_Not_Found, i;
             /* adding delvto */
 
 #ifdef RELAN
-            if (model->BSIM4type == -1)
+            if (model->BSIM4type == PMOS)
             {
-                if (here->BSIM4reliability->IsON != -1)
+                if (here->relStruct->IsON != -1)
                 {
-                    here->BSIM4delvto = here->BSIM4reliability->deltaVth ;
+                    here->BSIM4delvto = here->relStruct->deltaVth ;
                 }
             }
 #endif
diff --git a/src/spicelib/devices/bsim4/bsim4def.h b/src/spicelib/devices/bsim4/bsim4def.h
index 2d6b66252..f18846a4e 100644
--- a/src/spicelib/devices/bsim4/bsim4def.h
+++ b/src/spicelib/devices/bsim4/bsim4def.h
@@ -30,12 +30,7 @@ Modified by Pankaj Kumar Thakur, 07/23/2012
 #include "ngspice/noisedef.h"
 
 #ifdef RELAN
-typedef struct sBSIM4relStruct {
-    double time ;
-    double deltaVth ;
-    double t_star ;
-    int IsON ;
-} BSIM4relStruct ;
+#include "../relmodel/relmodeldefs.h"
 #endif
 
 typedef struct sBSIM4instance
@@ -584,7 +579,7 @@ typedef struct sBSIM4instance
 #endif /* NONOISE */
 
 #ifdef RELAN
-    BSIM4relStruct *BSIM4reliability ;
+    RELMODELrelStruct *relStruct ;
 #endif
 
 } BSIM4instance ;
@@ -837,6 +832,10 @@ typedef struct sBSIM4model
     BSIM4instance *BSIM4instances;
     IFuid BSIM4modName;
 
+#ifdef RELAN
+    GENmodel *BSIM4relmodelModel ;
+#endif
+
     /* --- end of generic struct GENmodel --- */
 
     int BSIM4type;
diff --git a/src/spicelib/devices/dev.c b/src/spicelib/devices/dev.c
index a4c7bc6b7..215090cc9 100644
--- a/src/spicelib/devices/dev.c
+++ b/src/spicelib/devices/dev.c
@@ -138,6 +138,10 @@ int add_udn(int,Evt_Udn_Info_t **);
 #include "ndev/ndevitf.h"
 #endif
 
+#ifdef RELAN
+#include "relmodel/relmodelitf.h"
+#endif
+
 static SPICEdev *(*static_devices[])(void) = {
     /* URC device MUST precede both resistors and capacitors */
     get_urc_info,
@@ -211,6 +215,10 @@ static SPICEdev *(*static_devices[])(void) = {
     get_ndev_info,
 #endif
 
+#ifdef RELAN
+    get_relmodel_info
+#endif
+
 };
 
 static int DEVNUM = NUMELEMS(static_devices);
diff --git a/src/spicelib/devices/mos1/Makefile.am b/src/spicelib/devices/mos1/Makefile.am
index c16190272..e523c113c 100644
--- a/src/spicelib/devices/mos1/Makefile.am
+++ b/src/spicelib/devices/mos1/Makefile.am
@@ -34,11 +34,6 @@ libmos1_la_SOURCES =	\
 	mos1trun.c
 
 
-if RELAN_WANTED
-libmos1_la_SOURCES += \
-	mos1reliability.c
-endif
-
 
 AM_CPPFLAGS = @AM_CPPFLAGS@ -I$(top_srcdir)/src/include
 AM_CFLAGS = $(STATIC)
diff --git a/src/spicelib/devices/mos1/mos1defs.h b/src/spicelib/devices/mos1/mos1defs.h
index 099006c04..c51cc2eb9 100644
--- a/src/spicelib/devices/mos1/mos1defs.h
+++ b/src/spicelib/devices/mos1/mos1defs.h
@@ -13,15 +13,6 @@ Modified: 2000 AlansFixes
 #include "ngspice/complex.h"
 #include "ngspice/noisedef.h"
 
-#ifdef RELAN
-typedef struct sMOS1relStruct {
-    double time ;
-    double deltaVth ;
-    double t_star ;
-    int IsON ;
-} MOS1relStruct ;
-#endif
-
 /* declarations for level 1 MOSFETs */
 
 /* information needed for each instance */
@@ -273,10 +264,6 @@ typedef struct sMOS1instance {
 #define MOS1dphibs_dw MOS1sens + 68
 #define MOS1dphibd_dw MOS1sens + 69
 
-#ifdef RELAN
-    MOS1relStruct *MOS1reliability ;
-#endif
-
 } MOS1instance ;
 
 #define MOS1vbd MOS1states+ 0   /* bulk-drain voltage */
diff --git a/src/spicelib/devices/mos1/mos1ext.h b/src/spicelib/devices/mos1/mos1ext.h
index dafe2f3bd..b1ed69d52 100644
--- a/src/spicelib/devices/mos1/mos1ext.h
+++ b/src/spicelib/devices/mos1/mos1ext.h
@@ -28,7 +28,3 @@ extern int MOS1convTest(GENmodel*,CKTcircuit*);
 extern int MOS1disto(int,GENmodel*,CKTcircuit*);
 extern int MOS1noise(int,int,GENmodel*,CKTcircuit*,Ndata*,double*);
 extern int MOS1dSetup(GENmodel*,CKTcircuit*);
-
-#ifdef RELAN
-extern int MOS1reliability (GENmodel *, CKTcircuit *, unsigned int) ;
-#endif
diff --git a/src/spicelib/devices/mos1/mos1init.c b/src/spicelib/devices/mos1/mos1init.c
index 6dc3f66db..a537f411c 100644
--- a/src/spicelib/devices/mos1/mos1init.c
+++ b/src/spicelib/devices/mos1/mos1init.c
@@ -73,12 +73,7 @@ SPICEdev MOS1info = {
  /* DEVacct       */ NULL,
 #endif                       
  /* DEVinstSize   */ &MOS1iSize,
- /* DEVmodSize    */ &MOS1mSize,
-
-#ifdef RELAN
- MOS1reliability
-#endif
-
+ /* DEVmodSize    */ &MOS1mSize
 };
 
 
diff --git a/src/spicelib/devices/mos1/mos1reliability.c b/src/spicelib/devices/mos1/mos1reliability.c
deleted file mode 100644
index 103d6b055..000000000
--- a/src/spicelib/devices/mos1/mos1reliability.c
+++ /dev/null
@@ -1,185 +0,0 @@
-/**********
-Author: 2015 Francesco Lannutti
-**********/
-
-#include "ngspice/ngspice.h"
-#include "ngspice/devdefs.h"
-#include "mos1defs.h"
-#include "ngspice/sperror.h"
-
-static int
-calculate_aging
-(
-    MOS1instance *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->MOS1reliability->t_star) / tau_0), beta)) * 1e9 <= 2)
-        {
-            here->MOS1reliability->deltaVth = Nt * A * pow (log (1 + pow (((t_aging + here->MOS1reliability->t_star) / tau_0), beta)), 2) ;
-        } else {
-            here->MOS1reliability->deltaVth = pow ((q / (4 * eps_0 * eps_hk)) * Nt * T_hk, 2) ;
-        }
-    } else {
-        here->MOS1reliability->deltaVth = here->MOS1reliability->deltaVth * log (1 + (1.718 / (1 + pow ((t_aging / tau_e), beta1)))) ;
-    }
-
-    if (!stress_or_recovery)
-    {
-        here->MOS1reliability->t_star = pow ((exp (sqrt (here->MOS1reliability->deltaVth / (Nt * A))) - 1), (1 / beta)) * tau_0 ;
-    }
-
-    return 0 ;
-}
-
-int
-MOS1reliability (GENmodel *inModel, CKTcircuit *ckt, unsigned int mode)
-{
-    MOS1model *model = (MOS1model *)inModel ;
-    MOS1instance *here ;
-    double delta, vds, vgs, von ;
-    int NowIsON ;
-
-    /*  loop through all the MOS1 device models */
-    for ( ; model != NULL ; model = model->MOS1nextModel)
-    {
-        /* loop through all the instances of the model */
-        for (here = model->MOS1instances ; here != NULL ; here=here->MOS1nextInstance)
-        {
-            if (model->MOS1type == -1)
-            {
-                vds = ckt->CKTstate0 [here->MOS1vds] ;
-                vgs = ckt->CKTstate0 [here->MOS1vgs] ;
-                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) ;
-                    if (vgs > von)
-                    {
-//                        printf ("Acceso!\n") ;
-                        NowIsON = 1 ;
-                    } else {
-//                        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) ;
-                    if (vgd > von)
-                    {
-//                        printf ("Acceso!\n") ;
-                        NowIsON = 1 ;
-                    } else {
-//                        printf ("Spento!\n") ;
-                        NowIsON = 0 ;
-                    }
-                }
-
-                // If it's the first time, initialize 'here->MOS1reliability->IsON'
-                if (here->MOS1reliability->IsON == -1)
-                {
-                    here->MOS1reliability->IsON = NowIsON ;
-                }
-
-                if (mode == 0)
-                {
-                    if (NowIsON)
-                    {
-                        if (here->MOS1reliability->IsON == 1)
-                        {
-                            // Until now, the device was ON - Do NOTHING
-                            delta = -1 ;
-                        } else if (here->MOS1reliability->IsON == 0) {
-                            // Until now, the device was OFF - Calculate recovery
-                            delta = ckt->CKTtime - here->MOS1reliability->time ;
-
-                            // Calculate Aging - Giogio Liatis' Model
-                            calculate_aging (here, delta, 0) ;
-
-                            // Update time and flag - Stress begins
-                            here->MOS1reliability->time = ckt->CKTtime ;
-                            here->MOS1reliability->IsON = 1 ;
-                        } else {
-                            fprintf (stderr, "Reliability Analysis Error\n") ;
-                        }
-                    } else {
-                        if (here->MOS1reliability->IsON == 1)
-                        {
-                            // Until now, the device was ON - Calculate stress
-                            delta = ckt->CKTtime - here->MOS1reliability->time ;
-
-                            // Calculate Aging - Giorgio Liatis' Model
-                            calculate_aging (here, delta, 1) ;
-
-                            // Update time and flag - Recovery begins
-                            here->MOS1reliability->time = ckt->CKTtime ;
-                            here->MOS1reliability->IsON = 0 ;
-                        } else if (here->MOS1reliability->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->MOS1reliability->IsON == 1)
-                    {
-                        // Calculate stress
-                        delta = ckt->CKTtime - here->MOS1reliability->time ;
-                        calculate_aging (here, delta, 1) ;
-
-                        // Update time and flag - Maybe Optional
-                        here->MOS1reliability->time = ckt->CKTtime ;
-                        here->MOS1reliability->IsON = 1 ;
-                    } else if (here->MOS1reliability->IsON == 0) {
-                        // Calculate recovery
-                        delta = ckt->CKTtime - here->MOS1reliability->time ;
-                        calculate_aging (here, delta, 0) ;
-
-                        // Update time and flag - Maybe Optional
-                        here->MOS1reliability->time = ckt->CKTtime ;
-                        here->MOS1reliability->IsON = 0 ;
-                    } 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") ;
-                }
-            }
-        }
-    }
-
-    return (OK) ;
-}
diff --git a/src/spicelib/devices/mos1/mos1set.c b/src/spicelib/devices/mos1/mos1set.c
index 334610cb4..521ff3fbf 100644
--- a/src/spicelib/devices/mos1/mos1set.c
+++ b/src/spicelib/devices/mos1/mos1set.c
@@ -206,17 +206,6 @@ do { if((here->ptr = SMPmakeElt(matrix, here->first, here->second)) == NULL){\
             TSTALLOC(MOS1SPbPtr,MOS1sNodePrime,MOS1bNode);
             TSTALLOC(MOS1SPdpPtr,MOS1sNodePrime,MOS1dNodePrime);
 
-#ifdef RELAN
-            if (model->MOS1type == -1)
-            {
-                here->MOS1reliability = TMALLOC (MOS1relStruct, 1) ;
-                here->MOS1reliability->time = 0 ;
-                here->MOS1reliability->deltaVth = 0 ;
-                here->MOS1reliability->t_star = 0 ;
-                here->MOS1reliability->IsON = -1 ;
-            }
-#endif
-
         }
     }
     return(OK);
diff --git a/src/spicelib/devices/relmodel/Makefile.am b/src/spicelib/devices/relmodel/Makefile.am
new file mode 100644
index 000000000..3cbdad580
--- /dev/null
+++ b/src/spicelib/devices/relmodel/Makefile.am
@@ -0,0 +1,23 @@
+## Process this file with automake to produce Makefile.in
+
+noinst_LTLIBRARIES = librelmodel.la
+
+if RELAN_WANTED
+librelmodel_la_SOURCES =	\
+	relmodel.c		\
+	relmodelcalcaging.c	\
+	relmodelmask.c		\
+	relmodelmpar.c		\
+	relmodeldefs.h		\
+	relmodelext.h		\
+	relmodelinit.c		\
+	relmodelinit.h		\
+	relmodelitf.h		\
+	relmodelsetup.c
+endif
+
+
+AM_CPPFLAGS = @AM_CPPFLAGS@ -I$(top_srcdir)/src/include
+AM_CFLAGS = $(STATIC)
+
+MAINTAINERCLEANFILES = Makefile.in
diff --git a/src/spicelib/devices/relmodel/relmodel.c b/src/spicelib/devices/relmodel/relmodel.c
new file mode 100644
index 000000000..9cafead6f
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodel.c
@@ -0,0 +1,35 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/devdefs.h"
+#include "relmodeldefs.h"
+#include "ngspice/suffix.h"
+
+IFparm RELMODELpTable [] = { /* parameters */
+} ;
+
+IFparm RELMODELmPTable [] = { /* model parameters */
+IOP ("k_b", RELMODEL_MOD_KB, IF_REAL, "Boltzmann Constant"),
+IOP ("h_cut", RELMODEL_MOD_HCUT, IF_REAL, "h_cut"),
+IOP ("nts", RELMODEL_MOD_NTS, IF_REAL, "Nts"),
+IOP ("eps_hk", RELMODEL_MOD_EPSHK, IF_REAL, "Dielectric Constant High-K"),
+IOP ("eps_SiO2", RELMODEL_MOD_EPSSIO2, IF_REAL, "Dielectric Constant SiO2"),
+IOP ("m_star", RELMODEL_MOD_MSTAR, IF_REAL, "m_star"),
+IOP ("w", RELMODEL_MOD_W, IF_REAL, "W"),
+IOP ("tau_0", RELMODEL_MOD_TAU0, IF_REAL, "tau_0"),
+IOP ("beta", RELMODEL_MOD_BETA, IF_REAL, "beta"),
+IOP ("tau_e", RELMODEL_MOD_TAUE, IF_REAL, "tau_e"),
+IOP ("beta1", RELMODEL_MOD_BETA1, IF_REAL, "beta1"),
+IP  ("relmodel", RELMODEL_MOD_RELMODEL,  IF_FLAG, "Flag to indicate RELMODEL")
+} ;
+
+char *RELMODELnames [] = {
+} ;
+
+int RELMODELnSize = NUMELEMS(RELMODELnames) ;
+int RELMODELpTSize = NUMELEMS(RELMODELpTable) ;
+int RELMODELmPTSize = NUMELEMS(RELMODELmPTable) ;
+int RELMODELiSize = 0 ;
+int RELMODELmSize = sizeof(RELMODELmodel) ;
diff --git a/src/spicelib/devices/relmodel/relmodelcalcaging.c b/src/spicelib/devices/relmodel/relmodelcalcaging.c
new file mode 100644
index 000000000..8d30c515c
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelcalcaging.c
@@ -0,0 +1,59 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/smpdefs.h"
+#include "ngspice/cktdefs.h"
+#include "relmodeldefs.h"
+#include "../bsim4/bsim4def.h"
+#include "ngspice/sperror.h"
+#include "ngspice/suffix.h"
+
+//#define CONSTepsZero (8.854214871e-12)   /* epsilon zero F/m */
+//#define CONSTepsSiO2 (3.4531479969e-11)  /* epsilon SiO2 F/m */
+
+int
+RELMODELcalculateAging (GENinstance *inInstance, int modType, double t_aging, unsigned int stress_or_recovery)
+{
+    double A, Nt ;
+    BSIM4instance *here ;
+    RELMODELmodel *relmodel ;
+
+    if (modType == 10)
+    {
+        here = (BSIM4instance *)inInstance ;
+    } else {
+        printf ("Error: The Reliability Analysis isn't supported for this model: %d\n\n\n", modType) ;
+        return 1 ;
+    }
+
+    relmodel = (RELMODELmodel *)(here->BSIM4modPtr->BSIM4relmodelModel) ;
+    if (relmodel == NULL)
+    {
+        printf ("Error: The Reliability Model hasn't been declared or it hasn't been connected to the device model (%s) through the '.appendmodel' command\n\n\n", here->BSIM4modPtr->BSIM4modName) ;
+        return 1 ;
+    }
+
+    Nt = pow ((sqrt (relmodel->RELMODELnts)), 3) * 1e-21 ;
+    A = (CHARGE / (4 * CONSTepsZero * 1e-9 * relmodel->RELMODELeps_hk)) * pow ((relmodel->RELMODELh_cut / (2 * sqrt (2 * relmodel->RELMODELm_star * relmodel->RELMODELw)) * 1e9), 2) ;
+
+    if (stress_or_recovery)
+    {
+        if (relmodel->RELMODELh_cut / (2 * sqrt (2 * relmodel->RELMODELm_star * relmodel->RELMODELw)) * log (1 + pow (((t_aging + here->relStruct->t_star) / relmodel->RELMODELtau_0), relmodel->RELMODELbeta)) * 1e9 <= 2)
+        {
+            here->relStruct->deltaVth = Nt * A * pow (log (1 + pow (((t_aging + here->relStruct->t_star) / relmodel->RELMODELtau_0), relmodel->RELMODELbeta)), 2) ;
+        } else {
+            here->relStruct->deltaVth = pow ((CHARGE / (4 * CONSTepsZero * 1e-9 * relmodel->RELMODELeps_hk)) * Nt * here->BSIM4modPtr->BSIM4toxe * 1e9, 2) ;
+        }
+    } else {
+        here->relStruct->deltaVth = here->relStruct->deltaVth * log (1 + (1.718 / (1 + pow ((t_aging / relmodel->RELMODELtau_e), relmodel->RELMODELbeta1)))) ;
+    }
+
+    if (!stress_or_recovery)
+    {
+        here->relStruct->t_star = pow ((exp (sqrt (here->relStruct->deltaVth / (Nt * A))) - 1), (1 / relmodel->RELMODELbeta)) * relmodel->RELMODELtau_0 ;
+    }
+
+    return 0 ;
+}
diff --git a/src/spicelib/devices/relmodel/relmodeldefs.h b/src/spicelib/devices/relmodel/relmodeldefs.h
new file mode 100644
index 000000000..a32ca9f59
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodeldefs.h
@@ -0,0 +1,71 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#ifndef RELMODEL
+#define RELMODEL
+
+#include "ngspice/ifsim.h"
+#include "ngspice/gendefs.h"
+#include "ngspice/cktdefs.h"
+#include "ngspice/complex.h"
+#include "ngspice/noisedef.h"
+
+typedef struct sRELMODELrelStruct {
+    double time ;
+    double deltaVth ;
+    double t_star ;
+    int IsON ;
+} RELMODELrelStruct ;
+
+typedef struct sRELMODELmodel
+{
+    int RELMODELmodType ;
+    struct sRELMODELmodel *RELMODELnextModel ;
+    IFuid RELMODELmodName ;
+
+    double RELMODELk_b ;
+    double RELMODELh_cut ;
+    double RELMODELnts ;
+    double RELMODELeps_hk ;
+    double RELMODELeps_SiO2 ;
+    double RELMODELm_star ;
+    double RELMODELw ;
+    double RELMODELtau_0 ;
+    double RELMODELbeta ;
+    double RELMODELtau_e ;
+    double RELMODELbeta1 ;
+    int RELMODELrelmodel ;
+
+    unsigned RELMODELk_bGiven : 1 ;
+    unsigned RELMODELh_cutGiven : 1 ;
+    unsigned RELMODELntsGiven : 1 ;
+    unsigned RELMODELeps_hkGiven : 1 ;
+    unsigned RELMODELeps_SiO2Given : 1 ;
+    unsigned RELMODELm_starGiven : 1 ;
+    unsigned RELMODELwGiven : 1 ;
+    unsigned RELMODELtau_0Given : 1 ;
+    unsigned RELMODELbetaGiven : 1 ;
+    unsigned RELMODELtau_eGiven : 1 ;
+    unsigned RELMODELbeta1Given : 1 ;
+    unsigned RELMODELrelmodelGiven : 1 ;
+} RELMODELmodel ;
+
+
+/* Global parameters */
+#define RELMODEL_MOD_KB       1
+#define RELMODEL_MOD_HCUT     2
+#define RELMODEL_MOD_NTS      3
+#define RELMODEL_MOD_EPSHK    4
+#define RELMODEL_MOD_EPSSIO2  5
+#define RELMODEL_MOD_MSTAR    6
+#define RELMODEL_MOD_W        7
+#define RELMODEL_MOD_TAU0     8
+#define RELMODEL_MOD_BETA     9
+#define RELMODEL_MOD_TAUE     10
+#define RELMODEL_MOD_BETA1    11
+#define RELMODEL_MOD_RELMODEL 12
+
+
+#include "relmodelext.h"
+#endif /* RELMODEL */
diff --git a/src/spicelib/devices/relmodel/relmodelext.h b/src/spicelib/devices/relmodel/relmodelext.h
new file mode 100644
index 000000000..827a8f74e
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelext.h
@@ -0,0 +1,9 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+extern int RELMODELcalculateAging (GENinstance *, int, double, unsigned int) ;
+extern int RELMODELmAsk (CKTcircuit *, GENmodel *, int, IFvalue *) ;
+extern int RELMODELmDelete (GENmodel **, IFuid, GENmodel *) ;
+extern int RELMODELmParam (int, IFvalue *, GENmodel *) ;
+extern int RELMODELsetup (SMPmatrix *, GENmodel *, CKTcircuit *, int *) ;
diff --git a/src/spicelib/devices/relmodel/relmodelinit.c b/src/spicelib/devices/relmodel/relmodelinit.c
new file mode 100644
index 000000000..503f2d292
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelinit.c
@@ -0,0 +1,70 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/config.h"
+
+#include "ngspice/devdefs.h"
+
+#include "relmodelitf.h"
+#include "relmodelext.h"
+#include "relmodelinit.h"
+
+
+SPICEdev RELMODELinfo = {
+    {
+        "RELMODEL",
+        "MOSFET Reliability Analysis Addiction Model",
+
+        &RELMODELnSize,
+        &RELMODELnSize,
+        RELMODELnames,
+
+        &RELMODELpTSize,
+        RELMODELpTable,
+
+        &RELMODELmPTSize,
+        RELMODELmPTable,
+
+        DEV_DEFAULT
+    },
+
+    NULL,           /* DEVparam       */
+    RELMODELmParam, /* DEVmodParam    */
+    NULL,           /* DEVload        */
+    RELMODELsetup,  /* DEVsetup       */
+    NULL,           /* DEVunsetup     */
+    NULL,           /* DEVpzSetup     */
+    NULL,           /* DEVtemperature */
+    NULL,           /* DEVtrunc       */
+    NULL,           /* DEVfindBranch  */
+    NULL,           /* DEVacLoad      */
+    NULL,           /* DEVaccept      */
+    NULL,           /* DEVdestroy     */
+    NULL,           /* DEVmodDelete   */
+    NULL,           /* DEVdelete      */
+    NULL,           /* DEVsetic       */
+    NULL,           /* DEVask         */
+    RELMODELmAsk,   /* DEVmodAsk      */
+    NULL,           /* DEVpzLoad      */
+    NULL,           /* DEVconvTest    */
+    NULL,           /* DEVsenSetup    */
+    NULL,           /* DEVsenLoad     */
+    NULL,           /* DEVsenUpdate   */
+    NULL,           /* DEVsenAcLoad   */
+    NULL,           /* DEVsenPrint    */
+    NULL,           /* DEVsenTrunc    */
+    NULL,           /* DEVdisto       */
+    NULL,           /* DEVnoise       */
+    NULL,           /* DEVsoaCheck    */
+    &RELMODELiSize, /* DEVinstSize    */
+    &RELMODELmSize, /* DEVmodSize     */
+    NULL           /* DEVreliability */
+} ;
+
+
+SPICEdev *
+get_relmodel_info (void)
+{
+    return &RELMODELinfo ;
+}
diff --git a/src/spicelib/devices/relmodel/relmodelinit.h b/src/spicelib/devices/relmodel/relmodelinit.h
new file mode 100644
index 000000000..e34d578fd
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelinit.h
@@ -0,0 +1,17 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#ifndef _RELMODELINIT_H
+#define _RELMODELINIT_H
+
+extern IFparm RELMODELpTable [] ;
+extern IFparm RELMODELmPTable [] ;
+extern char *RELMODELnames [] ;
+extern int RELMODELpTSize ;
+extern int RELMODELmPTSize ;
+extern int RELMODELnSize ;
+extern int RELMODELiSize ;
+extern int RELMODELmSize ;
+
+#endif
diff --git a/src/spicelib/devices/relmodel/relmodelitf.h b/src/spicelib/devices/relmodel/relmodelitf.h
new file mode 100644
index 000000000..2a12dffcb
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelitf.h
@@ -0,0 +1,10 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#ifndef DEV_RELMODEL
+#define DEV_RELMODEL
+
+SPICEdev *get_relmodel_info (void) ;
+
+#endif
diff --git a/src/spicelib/devices/relmodel/relmodelmask.c b/src/spicelib/devices/relmodel/relmodelmask.c
new file mode 100644
index 000000000..cc122d0dd
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelmask.c
@@ -0,0 +1,73 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/ifsim.h"
+#include "ngspice/cktdefs.h"
+#include "ngspice/devdefs.h"
+#include "relmodeldefs.h"
+#include "ngspice/sperror.h"
+#include "ngspice/suffix.h"
+
+int
+RELMODELmAsk (CKTcircuit *ckt, GENmodel *inst, int which, IFvalue *value)
+{
+    RELMODELmodel *model = (RELMODELmodel *)inst ;
+
+    NG_IGNORE (ckt) ;
+
+    switch (which)
+    {
+        case RELMODEL_MOD_KB:
+            value->rValue = model->RELMODELk_b ;
+            return (OK) ;
+
+        case RELMODEL_MOD_HCUT:
+            value->rValue = model->RELMODELh_cut ;
+            return (OK) ;
+
+        case RELMODEL_MOD_NTS:
+            value->rValue = model->RELMODELnts ;
+            return (OK) ;
+
+        case RELMODEL_MOD_EPSHK:
+            value->rValue = model->RELMODELeps_hk ;
+            return (OK) ;
+
+        case RELMODEL_MOD_EPSSIO2:
+            value->rValue = model->RELMODELeps_SiO2 ;
+            return (OK) ;
+
+        case RELMODEL_MOD_MSTAR:
+            value->rValue = model->RELMODELm_star ;
+            return (OK) ;
+
+        case RELMODEL_MOD_W:
+            value->rValue = model->RELMODELw ;
+            return (OK) ;
+
+        case RELMODEL_MOD_TAU0:
+            value->rValue = model->RELMODELtau_0 ;
+            return (OK) ;
+
+        case RELMODEL_MOD_BETA:
+            value->rValue = model->RELMODELbeta ;
+            return (OK) ;
+
+        case RELMODEL_MOD_TAUE:
+            value->rValue = model->RELMODELtau_e ;
+            return (OK) ;
+
+        case RELMODEL_MOD_BETA1:
+            value->rValue = model->RELMODELbeta1 ;
+            return (OK) ;
+
+        case RELMODEL_MOD_RELMODEL:
+            value->iValue = model->RELMODELrelmodel ;
+            return (OK) ;
+
+        default:
+            return (E_BADPARM) ;
+    }
+}
diff --git a/src/spicelib/devices/relmodel/relmodelmpar.c b/src/spicelib/devices/relmodel/relmodelmpar.c
new file mode 100644
index 000000000..8ba4cf8b0
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelmpar.c
@@ -0,0 +1,84 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "relmodeldefs.h"
+#include "ngspice/ifsim.h"
+#include "ngspice/sperror.h"
+#include "ngspice/suffix.h"
+#include "ngspice/const.h"
+
+int
+RELMODELmParam (int param, IFvalue *value, GENmodel *inMod)
+{
+    RELMODELmodel *model = (RELMODELmodel*)inMod ;
+
+    switch (param)
+    {
+        case RELMODEL_MOD_KB:
+            model->RELMODELk_b = value->rValue ;
+            model->RELMODELk_bGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_HCUT:
+            model->RELMODELh_cut = value->rValue ;
+            model->RELMODELh_cutGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_NTS:
+            model->RELMODELnts = value->rValue ;
+            model->RELMODELntsGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_EPSHK:
+            model->RELMODELeps_hk = value->rValue ;
+            model->RELMODELeps_hkGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_EPSSIO2:
+            model->RELMODELeps_SiO2 = value->rValue ;
+            model->RELMODELeps_SiO2Given = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_MSTAR:
+            model->RELMODELm_star = value->rValue ;
+            model->RELMODELm_starGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_W:
+            model->RELMODELw = value->rValue ;
+            model->RELMODELwGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_TAU0:
+            model->RELMODELtau_0 = value->rValue ;
+            model->RELMODELtau_0Given = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_BETA:
+            model->RELMODELbeta = value->rValue ;
+            model->RELMODELbetaGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_TAUE:
+            model->RELMODELtau_e = value->rValue ;
+            model->RELMODELtau_eGiven = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_BETA1:
+            model->RELMODELbeta1 = value->rValue ;
+            model->RELMODELbeta1Given = TRUE ;
+            return (OK) ;
+
+        case RELMODEL_MOD_RELMODEL:
+            model->RELMODELrelmodel = 1 ;
+            model->RELMODELrelmodelGiven = TRUE ;
+            break ;
+
+        default:
+            return (E_BADPARM) ;
+    }
+
+    return (OK) ;
+}
diff --git a/src/spicelib/devices/relmodel/relmodelsetup.c b/src/spicelib/devices/relmodel/relmodelsetup.c
new file mode 100644
index 000000000..91fe82dc7
--- /dev/null
+++ b/src/spicelib/devices/relmodel/relmodelsetup.c
@@ -0,0 +1,80 @@
+/**********
+Author: Francesco Lannutti - July 2015
+**********/
+
+#include "ngspice/ngspice.h"
+#include "ngspice/smpdefs.h"
+#include "ngspice/cktdefs.h"
+#include "relmodeldefs.h"
+#include "ngspice/sperror.h"
+#include "ngspice/suffix.h"
+
+int
+RELMODELsetup (SMPmatrix *matrix, GENmodel *inModel, CKTcircuit *ckt, int *states)
+{
+    RELMODELmodel *model = (RELMODELmodel *)inModel ;
+    NG_IGNORE (matrix) ;
+    NG_IGNORE (ckt) ;
+    NG_IGNORE (states) ;
+
+    /*  loop through all the RELMODEL device models */
+    for ( ; model != NULL ; model = model->RELMODELnextModel)
+    {
+        if (!model->RELMODELk_bGiven)
+        {
+            model->RELMODELk_b = 8.6e-5 ;
+        }
+
+        if (!model->RELMODELh_cutGiven)
+        {
+            model->RELMODELh_cut = 1.05e-34 ;
+        }
+
+        if (!model->RELMODELntsGiven)
+        {
+            model->RELMODELnts = 2e13 ;
+        }
+
+        if (!model->RELMODELeps_hkGiven)
+        {
+            model->RELMODELeps_hk = 25 ;
+        }
+
+        if (!model->RELMODELeps_SiO2Given)
+        {
+            model->RELMODELeps_SiO2 = 3.9 ;
+        }
+
+        if (!model->RELMODELm_starGiven)
+        {
+            model->RELMODELm_star = 0.1 * 9.11e-31 ;
+        }
+
+        if (!model->RELMODELwGiven)
+        {
+            model->RELMODELw = 1.5 * 1.6e-19 ;
+        }
+
+        if (!model->RELMODELtau_0Given)
+        {
+            model->RELMODELtau_0 = 1e-11 ;
+        }
+
+        if (!model->RELMODELbetaGiven)
+        {
+            model->RELMODELbeta = 0.373 ;
+        }
+
+        if (!model->RELMODELtau_eGiven)
+        {
+            model->RELMODELtau_e = 0.85e-9 ;
+        }
+
+        if (!model->RELMODELbeta1Given)
+        {
+            model->RELMODELbeta1 = 0.112 ;
+        }
+    }
+
+    return (OK) ;
+}
diff --git a/src/spicelib/parser/inp2dot.c b/src/spicelib/parser/inp2dot.c
index b2e7e5a72..bff6b955c 100644
--- a/src/spicelib/parser/inp2dot.c
+++ b/src/spicelib/parser/inp2dot.c
@@ -705,6 +705,14 @@ INP2dot(CKTcircuit *ckt, INPtables *tab, card *current, TSKtask *task, CKTnode *
 		/* don't have to do anything, since models were all done in
 		 * pass 1 */
 		goto quit;
+
+#ifdef RELAN
+        } else if (strcmp (token, ".appendmodel") == 0) {
+                /* don't have to do anything, since this command is handled in
+		 * pass 3 */
+                goto quit ;
+#endif
+
 	} else if ((strcmp(token, ".width") == 0) ||
 	           strcmp(token, ".print") == 0 || strcmp(token, ".plot") == 0) {
 		/* obsolete - ignore */
diff --git a/src/spicelib/parser/inp2m.c b/src/spicelib/parser/inp2m.c
index 038741380..90329be30 100644
--- a/src/spicelib/parser/inp2m.c
+++ b/src/spicelib/parser/inp2m.c
@@ -303,4 +303,12 @@ INP2M (CKTcircuit *ckt, INPtables * tab, card * current)
     if (waslead) {
         LITERR (" error:  no unlabeled parameter permitted on mosfet\n");
     }
+
+#ifdef RELAN
+    if (type == INPtypelook ("BSIM4"))
+    {
+        thismodel->INPmodfast->GENrelmodelModel = NULL ; /* Default value */
+    }
+#endif
+
 }
diff --git a/src/spicelib/parser/inpdomod.c b/src/spicelib/parser/inpdomod.c
index 954d2253b..e24cabb17 100644
--- a/src/spicelib/parser/inpdomod.c
+++ b/src/spicelib/parser/inpdomod.c
@@ -215,18 +215,39 @@ char *INPdomodel(CKTcircuit *ckt, card * image, INPtables * tab)
     /*  --------  Check if model is a MOSFET --------- */
     else if ((strcmp(type_name, "nmos") == 0)
 	       || (strcmp(type_name, "pmos") == 0)
+
+#ifdef RELAN
+	       || (strcmp (type_name, "relmodel") == 0)
+#endif
+
 	       || (strcmp(type_name, "nsoi") == 0)
 	       || (strcmp(type_name, "psoi") == 0)) {
 			err = INPfindLev(line, &lev);
 			switch (lev) {
 			case 0:
 			case 1:
-			    type = INPtypelook("Mos1");
-			    if (type < 0) {
-				err =
-				    INPmkTemp
-				    ("Device type MOS1 not available in this binary\n");
-			    }
+
+#ifdef RELAN
+                            if (!(strcmp (type_name, "relmodel")))
+                            {
+                                type = INPtypelook ("RELMODEL") ;
+                                if (type < 0)
+                                {
+                                    err = INPmkTemp ("Model type 'relmodel' not available in this binary\n") ;
+                                }
+                            } else {
+#endif
+
+                                type = INPtypelook ("Mos1") ;
+                                if (type < 0)
+                                {
+                                    err = INPmkTemp ("Device type MOS1 not available in this binary\n") ;
+                                }
+
+#ifdef RELAN
+                            }
+#endif
+
 			    break;
 			case 2:
 			    type = INPtypelook("Mos2");
diff --git a/src/spicelib/parser/inppas3.c b/src/spicelib/parser/inppas3.c
index 7fdbda01e..8f0fa8920 100644
--- a/src/spicelib/parser/inppas3.c
+++ b/src/spicelib/parser/inppas3.c
@@ -14,6 +14,121 @@ Modified: AlansFixes
 
 extern IFsimulator *ft_sim;
 
+#ifdef RELAN
+static void
+dot_appendmodel (CKTcircuit *ckt, INPtables *tab, card *current)
+{
+    /* .appendmodel from_model to_model */
+
+    int type;                  /* the type the model says it is */
+    char *line ;               /* the part of the current line left to parse */
+    char *relmodel_name ;      /* the relmodel name */
+    char *appendmodel_name ;   /* the appendmodel name */
+    char *save ;               /* saj - used to save the posn of the start of
+                                  the parameters if the model is a mosfet */
+    char *model ;              /* the name of the model */
+    INPmodel *thismodel ;      /* pointer to model description for user's model */
+    GENmodel *mdfast ;         /* pointer to the actual model */
+    char* err_msg ;
+
+#ifdef TRACE
+    printf("INP2M: Parsing '%s'\n",current->line);
+#endif
+
+    line = current->line ;
+
+    INPgetTok (&line, &relmodel_name, 1) ; /* throw away .appendmodel */
+    tfree (relmodel_name) ;
+
+    INPgetTok (&line, &relmodel_name, 1) ; /* get relmodel name */
+
+    save = line ; /* saj - save the posn for later if the default mosfet model is used */
+    thismodel = NULL ;
+
+    err_msg = INPgetMod (ckt, relmodel_name, &thismodel, tab) ;
+    tfree (err_msg) ;
+
+    /* check if using model binning -- pass in line since need 'l' and 'w' */
+    if (thismodel == NULL)
+    {
+        INPgetModBin (ckt, relmodel_name, &thismodel, tab, line) ;
+    }
+
+    model = relmodel_name ; /* mosfet */
+    line = save ;        /* reset the posn to what it sould be */
+
+    INPinsert (&model, tab) ;
+    thismodel = NULL ;
+
+#ifdef TRACE
+    printf("INP2M: Looking up model\n");
+#endif
+
+    err_msg = INPgetMod (ckt, model, &thismodel, tab) ;
+    if (thismodel == NULL)
+    {
+        INPgetModBin (ckt, model, &thismodel, tab, save) ;
+        if (thismodel == NULL)
+        {
+            current->error = err_msg ;
+        } else {
+            tfree (err_msg) ;
+        }
+    } else {
+        tfree (err_msg) ;
+    }
+
+    if (thismodel != NULL)
+    {
+        if (thismodel->INPmodType != INPtypelook ("Mos1")
+                && thismodel->INPmodType != INPtypelook ("Mos2")
+                && thismodel->INPmodType != INPtypelook ("Mos3")
+                && thismodel->INPmodType != INPtypelook ("Mos5")
+                && thismodel->INPmodType != INPtypelook ("Mos6")
+                && thismodel->INPmodType != INPtypelook ("Mos8")
+                && thismodel->INPmodType != INPtypelook ("Mos9")
+                && thismodel->INPmodType != INPtypelook ("BSIM1")
+                && thismodel->INPmodType != INPtypelook ("BSIM2")
+                && thismodel->INPmodType != INPtypelook ("BSIM3")
+                && thismodel->INPmodType != INPtypelook ("BSIM3v32")
+                && thismodel->INPmodType != INPtypelook ("B4SOI")
+                && thismodel->INPmodType != INPtypelook ("B3SOIPD")
+                && thismodel->INPmodType != INPtypelook ("B3SOIFD")
+                && thismodel->INPmodType != INPtypelook ("B3SOIDD")
+                && thismodel->INPmodType != INPtypelook ("BSIM4")
+                && thismodel->INPmodType != INPtypelook ("BSIM4v5")
+                && thismodel->INPmodType != INPtypelook ("BSIM4v6")
+                && thismodel->INPmodType != INPtypelook ("BSIM4v7")
+                && thismodel->INPmodType != INPtypelook ("BSIM3v0")
+                && thismodel->INPmodType != INPtypelook ("BSIM3v1")
+                && thismodel->INPmodType != INPtypelook ("SOI3")
+                && thismodel->INPmodType != INPtypelook ("HiSIM2")
+                && thismodel->INPmodType != INPtypelook ("HiSIMHV")
+                && thismodel->INPmodType != INPtypelook ("RELMODEL")
+           )
+        {
+            LITERR ("incorrect model type") ;
+            return ;
+        }
+        type = thismodel->INPmodType ;
+        mdfast = (thismodel->INPmodfast) ;
+
+        INPgetTok (&line, &appendmodel_name, 1) ; /* get appendmodel name */
+        thismodel = NULL ;
+        err_msg = INPgetMod (ckt, appendmodel_name, &thismodel, tab) ;
+        tfree (err_msg) ;
+        if (thismodel == NULL)
+        {
+            INPgetModBin (ckt, relmodel_name, &thismodel, tab, line) ;
+        }
+        thismodel->INPmodfast->GENrelmodelModel = mdfast ; /* Attach 'relmodel' model to 'appendmodel' model */
+    } else {
+        LITERR ("Error: The model 'relmodel' hasn't been defined\n\n") ;
+    }
+
+    return ;
+}
+#endif
 
 /* pass 3 - Read all nodeset and IC lines. All circuit nodes will have
  * been created by now, (except for internal device nodes), so any
@@ -144,6 +259,14 @@ INPpas3(CKTcircuit *ckt, card *data, INPtables *tab, TSKtask *task,
                 break;
             }
         }
+
+#ifdef RELAN
+        else if (strcmp (token, ".appendmodel") == 0) {
+                dot_appendmodel (ckt, tab, current) ;
+                goto quit ;
+        }
+#endif
+
     }
 quit:
     FREE(token);