only zero the Jacobian once (#1636)

now it is the role of the jac() wrapper to zero the analytic Jacobian.
The networks can assume that it is zero.  Previously both VODE
and the network were doing the zeroing, which was redundant.

integration/VODE/vode_dvjac.H

@@ -79,9 +79,6 @@ void dvjac (int& IERPJ, BurnT& state, DvodeT& vstate)
             // Indicate that the Jacobian is current for this solve.
             vstate.JCUR = 1;
 
-            // Initialize the Jacobian to zero
-            vstate.jac.zero();
-
             jac(vstate.tn, state, vstate, vstate.jac);
 
 #ifdef ALLOW_JACOBIAN_CACHING

sphinx_docs/source/integrators.rst

@@ -276,6 +276,11 @@ The form looks like:
    A network is not required to provide a Jacobian if a numerical
    Jacobian is used.
 
+.. important::
+
+   The integrator does not zero the Jacobian elements.  It is the responsibility
+   of the Jacobian implementation to zero the Jacobian array if necessary.
+
 
 Jacobian wrapper
 ^^^^^^^^^^^^^^^^