I am currently working on developing a compact model of a saturating FET IV in Verilog-A. All of the governing equations are correct and they have shown great results in MATLAB. In Cadence, when simulating a single n-type or p-type device, DC analysis works fine. When performing transient analysis, it's possible to simulate a NAND and an Inverter with correct results. These simulations tend to take a very long time, possibly due to multiple necessary iterations to achieve convergence. The problem is that I cannot seem to get the DC analysis to converge when simulating a simple two-device inverter (one n-type, one p-type). I have adjusted the simulation convergence parameters with no success. I replaced every "exp" function in the model with the "limexp" function with no success. We also tried to use the built in convergence functions (slew, transition, etc.). Now I strongly believe the problem is due to a discontinuity in the 1st derivative of the IV curve. I've tried to announce it in the model using the discontinuity function, but it did not seem to have any effect. It is possible that I did not implement this function properly, however.
In MATLAB I used convolution to smooth out the curve with great results. Is there a way to achieve this in Verilog-A? Is there any other possible ways to smoothen the output IV so that the simulator can converge? Do you have any other suggestions?

Here is the error code from the DC inverter simulation:

Circuit inventory:
nodes 3
capacitor 1
gnr 2
vsource 2

Entering remote command mode using MPSC service (spectre, ipi, v0.0, spectre0_9383_1, ).

*********************************
DC Analysis `dc': Vmax = (0 -> 1)
*********************************
Important parameter values:
reltol = 1e-03
abstol(I) = 1 pA
abstol(V) = 1 uV
temp = 27 C
tnom = 27 C
tempeffects = all
gmin = 1 pS
Trying `homotopy = gmin'.
Trying `homotopy = source'.
Trying `homotopy = gmin'.
Trying `homotopy = source'.
Trying `homotopy = dptran'.
DPTran step# = 1000
DPTran step# = 2000
DPTran step# = 3000
DPTran step# = 4000
DPTran step# = 5000
DPTran step# = 6000
DPTran step# = 7000
DPTran step# = 8000
DPTran step# = 9000
DPTran step# = 10000
DPTran step# = 11000
DPTran step# = 12000
DPTran step# = 13000
DPTran step# = 14000
DPTran step# = 15000
Step limit exceeded.
Trying `homotopy = ptran'.
Trying `homotopy = arclength'.
None of the instantiated devices support arclength homotopy. Skipping.

Arclength failed after a total of 0 iterations.
Attempting Newton loop with latest solution.

Error found by spectre at Vmax = 20e-03 during DC analysis `dc'.
ERROR (SPECTRE-16080): No DC solution found (no convergence). Last acceptable solution computed at 0.

The values for those nodes that did not converge on the last Newton iteration are given below. The manner in which the convergence criteria were not satisfied is also given.
Failed test: | Value | > RelTol*Ref + AbsTol

V(net016) = 0 V, previously 1 V.
residue too large: | 623.973 nA | > 622.808 nA + 1 pA
V(Out) = 0 V, previously 559.927 mV.
residue too large: | -623.973 nA | > 622.808 nA + 1 pA
I(V0:p) = 0 A, previously -1.64376 mA.
update too large: | -5.67325 pA | > 0 A + 1 pA


Thank you for your help

I am doing this work for a class so I have to work on the given equations with a discontinuity in the 2nd derivative at the edge of the sub-threshold. I believe all the compact model equations will have a some kind of discontinuity in their IVs since at that point equations change from the sub-threshold zone to linear zone. The slopes of these equations would be completely different. Somehow in complex Verilog-A model files like EKV model they are by passing this discontinuity problem. I have been reading about the discontinuity function, but I can not make it work. My individual n-type and p-type devices work fine, I have the convergence error only when I try to simulate an inverter with DC sweep.  Thanks for your help!