Countercurrent spontaneous imbibition (SI) experiments are among the most common multiphase experiments performed on porous rock samples. Although the samples rarely are designed to give linear flow, they are often modeled and interpreted using mathematical descriptions assuming linear flow. In this work, the goal is to derive general understanding of how imbibition into different sample geometries behaves compared with linear (1D) imbibition.

Using the mathematical theory of $N$-volume spheres ($N$ being the dimension), we consider core samples as $N$-spherical and quantify their geometry by a dimension $N$ and length scale $L$. For the special cases, $N=1,2,3$, we obtain linear, radial, and spherical flow; however, we treat the dimension as an arbitrary real number for cases not adhering to either of these regimes. Particularly, for rectangular or cylindrical core plugs, a continuous range of dimensions is produced. Theoretical calculations of $N$ and $L$ of arbitrary sample shapes are derived based on relations with area per volume and derivative of area with respect to volume. They correctly produce limit cases and physically meaningful values for symmetrical, cylindrical, and rectangular geometries. The differential equation for countercurrent imbibition in $N$-dimensions is derived and solved with numerical examples. Also, a simplified analytical solution assuming piston-like displacement is derived to get illustrative relations between system parameters (including dimension) and recovery and front position.

Predicted recovery profiles of the $N$-dimensional (N-D) solution overlap consistently with numerical simulations (by an alternative simulator) into cylindrical geometries with a wide range of height/diameter ratios and viscosity ratios. At early time, the saturation profiles are self-similar (look the same plotted against position divided by square root of time) and identical regardless of dimension. As long as the profiles are self-similar near the open boundary, recovery is proportional to the square root of time. For $N=1$, this lasts long after the front has reached the closed boundary, while for N >1, it can happen long before the boundary has been reached. The same time scale was applicable for all geometries, stating especially that doubling the length scale increased the time of the entire recovery profile by a factor 4. As long as recovery is proportional to square root of time, at a given time, it is also proportional to dimension over length, $NL$, and the time needed to reach same level of recovery is proportional to $LN2$. Permeability anisotropy (lower vertical than horizontal permeability) could effectively be modeled using an effective increased height, which further was captured by the dimension and length scale. Literature and in-house experimental data were matched by the model and used to validate model predictions such as the variation in time, shape of recovery curves with changes in dimension, and the importance of accounting for the dimension (geometry) during data interpretation. The model is valid for all wetting states but assumes negligible gravity and compressibility effects.