nmodl icon indicating copy to clipboard operation
nmodl copied to clipboard

MOD files containing two SOLVEs don't compile.

Open 1uc opened this issue 1 year ago • 1 comments

For example the following:

$ cat two_blocks.mod
NEURON {
    SUFFIX two_blocks
    RANGE A, B, X, Y
    NONSPECIFIC_CURRENT il
}

STATE {
    X
    Y
    A
    B
}

ASSIGNED {
    il
    i
}

BREAKPOINT {
    SOLVE state1 METHOD sparse
    SOLVE state2 METHOD sparse
    il = i
}

INITIAL {
    A = 1.0
    B = 3.0
    X = 0.0
    Y = 1.0
}

KINETIC state1 {
    ~ X <-> Y (0.4, 0.5)
}

KINETIC state2 {
    ~ A <-> B (0.6, 0.7)
    i = (f_flux-b_flux)
}

generates the following C++ code:

            Eigen::Matrix<double, 2, 1> nmodl_eigen_xm;
            double* nmodl_eigen_x = nmodl_eigen_xm.data();
            nmodl_eigen_x[static_cast<int>(0)] = inst->X[id];
            nmodl_eigen_x[static_cast<int>(1)] = inst->Y[id];
            // call newton solver
            functor_two_blocks_0 newton_functor(nt, inst, id, pnodecount, v, indexes, data, thread);
            newton_functor.initialize();
            int newton_iterations = nmodl::newton::newton_solver(nmodl_eigen_xm, newton_functor);
            if (newton_iterations < 0) assert(false && "Newton solver did not converge!");
            inst->X[id] = nmodl_eigen_x[static_cast<int>(0)];
            inst->Y[id] = nmodl_eigen_x[static_cast<int>(1)];
            newton_functor.initialize(); // TODO mimic calling F again.
            newton_functor.finalize();

            
            Eigen::Matrix<double, 2, 1> nmodl_eigen_xm;
            double* nmodl_eigen_x = nmodl_eigen_xm.data();
            nmodl_eigen_x[static_cast<int>(0)] = inst->A[id];
            nmodl_eigen_x[static_cast<int>(1)] = inst->B[id];
            // call newton solver
            functor_two_blocks_1 newton_functor(nt, inst, id, pnodecount, v, indexes, data, thread);
            newton_functor.initialize();
            int newton_iterations = nmodl::newton::newton_solver(nmodl_eigen_xm, newton_functor);
            if (newton_iterations < 0) assert(false && "Newton solver did not converge!");
            inst->A[id] = nmodl_eigen_x[static_cast<int>(0)];
            inst->B[id] = nmodl_eigen_x[static_cast<int>(1)];
            newton_functor.initialize(); // TODO mimic calling F again.
            newton_functor.finalize();

and fails with:

$ nrnivmodl -coreneuron
x86_64/corenrn/mod2c/two_blocks.cpp:786:41: error: redeclaration of ‘Eigen::Matrix<double, 2, 1> nmodl_eigen_xm’
  786 |             Eigen::Matrix<double, 2, 1> nmodl_eigen_xm;
      |                                         ^~~~~~~~~~~~~~
x86_64/corenrn/mod2c/two_blocks.cpp:771:41: note: ‘Eigen::Matrix<double, 2, 1> nmodl_eigen_xm’ previously declared here
  771 |             Eigen::Matrix<double, 2, 1> nmodl_eigen_xm;
      |                                         ^~~~~~~~~~~~~~

NOCMODL seems to generate code that compiles and looks sane. Therefore, it's not obvious the MOD file is invalid.

1uc avatar Oct 07 '24 11:10 1uc

Summary of discussion with @nrnhines:

  1. Carefully check if NOCMODL supports this and check ModelDB CI for usage examples.
  2. Also consider recursive use, e.g. solving a LINEAR block inside a KINETIC block.
  3. Note that STATE variables are declared globally, i.e. there's no notion of the STATEs of block state1. Some solvers might require knowing the number of states.

Generally, this seems to be a problem that can be revisited when there's an instance of this problem.

1uc avatar Oct 07 '24 12:10 1uc