fix: ww crbl demo updates

Author: maedoc

tvb_documentation/demos/ww-crbl/run_ww_crbl_hybrid.py

@@ -209,6 +209,7 @@ def _get_supported_models_classes() -> tuple:
         ReducedSetFitzHughNagumo,
         ReducedSetHindmarshRose,
     )
+    from tvb.simulator.models.cerebellar_mf import CerebellarMF
     from tvb.simulator.models.linear import Linear
     _SUPPORTED_MODELS_CACHE = (
         MontbrioPazoRoxin,
@@ -236,6 +237,7 @@ def _get_supported_models_classes() -> tuple:
         ZetterbergJansen,
         ReducedSetFitzHughNagumo,
         ReducedSetHindmarshRose,
+        CerebellarMF,
         Linear,
     )
     return _SUPPORTED_MODELS_CACHE
@@ -1486,8 +1488,14 @@ def _run_compiled(
                 args.append(stim_arr)
 
         # Per-subnetwork spatial parameter arrays (heterogeneous per-node parameters)
+        # Custom-template models bake all params into generated code, so skip
+        # spatial_parameter_names auto-detection (which may pick up non-spatial
+        # arrays like P_grc[5] that can't broadcast to n_nodes).
         for sn_info in analysis.subnetworks:
-            sp_names = list(getattr(sn_info.model, 'spatial_parameter_names', []))
+            if hasattr(sn_info.model, '_nb_hybrid_custom_template'):
+                sp_names = []
+            else:
+                sp_names = list(getattr(sn_info.model, 'spatial_parameter_names', []))
             if sp_names:
                 sp_arr = np.array(
                     [np.broadcast_to(getattr(sn_info.model, n), (sn_info.n_nodes,)).ravel()

tvb_documentation/demos/ww-crbl/ww_crbl_fc.png

@@ -295,13 +295,15 @@ def compute_coupling_${p.name}(
                "heun" if isinstance(sn.integrator, HeunDeterministic) else "euler"
     svars = list(sn.model.state_variables)
 
-    # Zerlaut models use a custom template for dfun generation;
-    # they don't have coupling_terms / state_variable_dfuns / global_parameter_names.
+    # Models with custom Mako templates (ZerlautAdaptation*, CerebellarMF)
+    # don't have state_variable_dfuns / global_parameter_names, but they do
+    # have coupling_terms.
     _has_custom_template = hasattr(sn.model, '_nb_hybrid_custom_template')
     if _has_custom_template:
-        cterms = ['Coupling_Term']  # Zerlaut: single cvar=[0]
+        cterms = list(sn.model.coupling_terms)
         dfuns = None
         gparams = {}
+        sparams_list = []  # custom template bakes all params — no spatial params needed
     else:
         cterms = list(sn.model.coupling_terms)
         dfuns = sn.model.state_variable_dfuns
@@ -312,16 +314,16 @@ def compute_coupling_${p.name}(
     n_nodes = sn.n_nodes
     n_modes = sn.n_modes
     svb = sn.model.state_variable_boundaries
-    lo_map = {k: float(v[0]) if v[0] != float('-inf') and not (v[0] != v[0]) else None for k, v in svb.items()} if svb else {}
-    hi_map = {k: float(v[1]) if v[1] != float('inf') and not (v[1] != v[1]) else None for k, v in svb.items()} if svb else {}
+    lo_map = {}
+    hi_map = {}
     import math as _math, numpy as _np
     lo_map = {}
     hi_map = {}
     if svb:
         import numpy as _np
         for k, v in svb.items():
-            lo_map[k] = float(v[0]) if _np.isfinite(v[0]) else None
-            hi_map[k] = float(v[1]) if _np.isfinite(v[1]) else None
+            lo_map[k] = float(v[0]) if v[0] is not None and _np.isfinite(v[0]) else None
+            hi_map[k] = float(v[1]) if v[1] is not None and _np.isfinite(v[1]) else None
     svars_str = ', '.join(svars)
     cterms_str = ', '.join(cterms)
     i1svars_str = ', '.join(['i1' + s for s in svars])

tvb_documentation/demos/ww-crbl/ww_crbl_hybrid_demo.ipynb

@@ -428,9 +428,11 @@ class CerebellarMF(ModelNumbaDfun):
 
     external_input_ex_ex = NArray(
         label=":math:`\\nu_e^{drive}`",
-        default=numpy.array([0.315e-3]),
-        domain=Range(lo=0.0, hi=0.1, step=0.001),
-        doc="External excitatory drive to excitatory population (GrC/GoC)")
+        default=numpy.array([0.05]),
+        domain=Range(lo=0.0, hi=0.5, step=0.001),
+        doc="External excitatory drive to excitatory population (GrC/GoC). "
+            "The default 0.05 produces oscillatory dynamics with peak/mean > 3x "
+            "in theta-alpha range across connected CRBL nodes.")
 
     external_input_ex_in = NArray(
         label=":math:`\\nu_i^{drive}`",
@@ -464,33 +466,39 @@ class CerebellarMF(ModelNumbaDfun):
 
     mf_to_grc = NArray(
         label="Mossy→GrC split",
-        default=numpy.array([0.97]),
+        default=numpy.array([1.0]),
         domain=Range(lo=0.0, hi=1.0, step=0.01),
-        doc="Fraction of mossy fiber input received by GrC")
+        doc="Fraction of mossy fiber input received by GrC. "
+            "Set to 1.0 to match the multimf_ww monolithic model where "
+            "c_mossy is the full coupling signal without anatomical sub-fractions.")
 
     mf_to_goc = NArray(
         label="Mossy→GoC split",
-        default=numpy.array([0.03]),
+        default=numpy.array([1.0]),
         domain=Range(lo=0.0, hi=1.0, step=0.01),
-        doc="Fraction of mossy fiber input received by GoC")
+        doc="Fraction of mossy fiber input received by GoC. "
+            "Set to 1.0 to match the multimf_ww monolithic model.")
 
     pf_to_goc = NArray(
         label="Parallel→GoC split",
-        default=numpy.array([0.14]),
+        default=numpy.array([1.0]),
         domain=Range(lo=0.0, hi=1.0, step=0.01),
-        doc="Fraction of parallel fiber input received by GoC")
+        doc="Fraction of parallel fiber input received by GoC. "
+            "Set to 1.0 to match the multimf_ww monolithic model.")
 
     pf_to_mli = NArray(
         label="Parallel→MLI split",
-        default=numpy.array([0.55]),
+        default=numpy.array([1.0]),
         domain=Range(lo=0.0, hi=1.0, step=0.01),
-        doc="Fraction of parallel fiber input received by MLI")
+        doc="Fraction of parallel fiber input received by MLI. "
+            "Set to 1.0 to match the multimf_ww monolithic model.")
 
     pf_to_pc = NArray(
         label="Parallel→PC split",
-        default=numpy.array([0.31]),
+        default=numpy.array([1.0]),
         domain=Range(lo=0.0, hi=1.0, step=0.01),
-        doc="Fraction of parallel fiber input received by PC")
+        doc="Fraction of parallel fiber input received by PC. "
+            "Set to 1.0 to match the multimf_ww monolithic model.")
 
     # -----------------------------------------------------------------------
     # Model metadata
@@ -499,6 +507,9 @@ class CerebellarMF(ModelNumbaDfun):
         label="Coupling terms",
         default=["mossy", "parallel"])
 
+    # nb-hybrid: use custom CerebellarMF Mako template — do NOT add state_variable_dfuns
+    _nb_hybrid_custom_template = "nb-cerebellar-dfun.py.mako"
+
     parameter_names = List(
         of=str,
         label="List of parameters for this model",
@@ -527,10 +538,10 @@ class CerebellarMF(ModelNumbaDfun):
     state_variable_range = Final(
         label="State Variable ranges [lo, hi]",
         default={
-            "GrC":   numpy.array([0.5e3, 0.5e3]),
-            "GoC":   numpy.array([5.0e3, 5.0e3]),
-            "MLI":   numpy.array([15.0e3, 15.0e3]),
-            "PC":    numpy.array([38.0e3, 38.0e3]),
+            "GrC":   numpy.array([0.1, 0.1]),
+            "GoC":   numpy.array([0.02, 0.02]),
+            "MLI":   numpy.array([0.2, 0.2]),
+            "PC":    numpy.array([0.1, 0.1]),
             "noise": numpy.array([0.0, 0.0]),
         },
         doc="""Expected dynamic range for each state variable.
@@ -666,7 +677,7 @@ def TF_inhibitory_goc(self, fe, fi, fe_ext, fi_ext, W=0):
         return self._TF_goc(fe, fi, fe_ext, fi_ext, W,
                             self.P_goc, self.Q_grc_goc, self.Q_goc_goc,
                             self.tau_grc_goc, self.tau_goc_goc,
-                            self.E_i, self.E_i,
+                            self.E_e, self.E_i,
                             self.g_L_goc, self.C_m_goc, self.E_L_goc,
                             self.K_grc_goc, self.K_goc_goc,
                             self.Q_mf_goc, self.tau_mf_goc,