change volume method (#2781)

share/lib/python/neuron/rxd/region.py

@@ -344,8 +344,14 @@ def __repr__(self):
     def _short_repr(self):
         return "Extracellular"
 
-    def volume(self, index):
+    def volume(self, index=None):
         """Returns the volume of the voxel at a given index"""
+        if index is None:
+            if numpy.isscalar(self.alpha):
+                vol = self._nx * self._ny * self._nz * numpy.prod(self._dx) * self.alpha
+            else:
+                vol = numpy.sum(self.alpha) * numpy.prod(self._dx)
+            return vol
         if numpy.isscalar(self.alpha):
             return numpy.prod(self._dx) * self.alpha
         return numpy.prod(self._dx) * self.alpha[index]
@@ -1000,6 +1006,16 @@ def secs(self, value):
         else:
             raise RxDException("Cannot set secs now; model already instantiated")
 
-    def volume(self, index):
+    def volume(self, index=None):
         """Returns the volume of the voxel at a given index"""
+        initializer._do_init()
+        if index is None:
+            vol = 0
+            if hasattr(self, "_vol") and any(self._secs3d):
+                vol += numpy.sum(self._vol)
+            if hasattr(self, "_geometry") and any(self._secs1d):
+                vol += numpy.sum(
+                    [self._geometry.volumes1d(sec) for sec in self._secs1d]
+                )
+            return vol
         return self._vol[index]

share/lib/python/neuron/rxd/rxd.py

@@ -905,25 +905,27 @@ def _setup_matrices():
             num_3d_indices_per_1d_seg = numpy.asarray(
                 num_3d_indices_per_1d_seg, dtype=numpy.int64
             )
+            hybrid_grid_ids.append(-1)
             hybrid_grid_ids = numpy.asarray(hybrid_grid_ids, dtype=numpy.int64)
 
             hybrid_indices3d = numpy.asarray(hybrid_indices3d, dtype=numpy.int64)
             rates = numpy.asarray(rates, dtype=float)
             volumes1d = numpy.asarray(volumes1d, dtype=float)
             volumes3d = numpy.asarray(volumes3d, dtype=float)
             dxs = numpy.asarray(grids_dx, dtype=float)
-            set_hybrid_data(
-                num_1d_indices_per_grid,
-                num_3d_indices_per_grid,
-                hybrid_indices1d,
-                hybrid_indices3d,
-                num_3d_indices_per_1d_seg,
-                hybrid_grid_ids,
-                rates,
-                volumes1d,
-                volumes3d,
-                dxs,
-            )
+            if hybrid_grid_ids.size > 1:
+                set_hybrid_data(
+                    num_1d_indices_per_grid,
+                    num_3d_indices_per_grid,
+                    hybrid_indices1d,
+                    hybrid_indices3d,
+                    num_3d_indices_per_1d_seg,
+                    hybrid_grid_ids,
+                    rates,
+                    volumes1d,
+                    volumes3d,
+                    dxs,
+                )
 
         # TODO: Replace this this to handle 1d/3d hybrid models
         """

test/rxd/test_region_volume.py

@@ -0,0 +1,20 @@
+def test_region_volume(neuron_nosave_instance):
+    h, rxd, _ = neuron_nosave_instance
+    dend1 = h.Section("dend1")
+    dend2 = h.Section("dend2")
+    diff = 1e-10
+    dend1.nseg = 4
+    dend2.nseg = 4
+    cyt1 = rxd.Region(dend1.wholetree(), "i", dx=0.25)
+    cyt2 = rxd.Region(dend2.wholetree(), "i", dx=0.25)
+    ca1 = rxd.Species(cyt1, name="ca1", charge=2, initial=1e-12)
+    ca2 = rxd.Species(cyt2, name="ca2", charge=2, initial=1e-12)
+    dend1.L = 10
+    dend1.diam = 2
+    dend2.L = 10
+    dend2.diam = 2
+    rxd.set_solve_type(domain=[dend1], dimension=1)
+    rxd.set_solve_type(domain=[dend2], dimension=3)
+    h.finitialize()
+    assert abs(cyt1.volume() - sum(ca1.nodes(cyt1).volume)) < diff
+    assert abs(cyt2.volume() - sum(ca2.nodes(cyt2).volume)) < diff