Multiphysics Modeling of a Hybrid Photovoltaic-Plant Microbial Fuel Cell (Biopile) System: Numerical Predictions versus Field Measurements ()
1. Introduction
This hybrid off-grid and mini-grid architectures in tropical, resource-constrained environments increasingly combine intermittent photovoltaics (PV) with complementary base-load sources to minimize battery storage requirements and secure service continuity. Plant microbial fuel cells (P-MFCs) offer a promising complementary pathway: rhizodeposits excreted by living plants continuously nourish electroactive microbial communities in anoxic soil environments, sustaining power generation through biological and mass-transport dynamics independent of instantaneous solar irradiance [1]-[3]. Building upon conventional microbial fuel cell (MFC) mechanics [4] [5], P-MFCs convert light energy into persistent bio-electricity via living root networks [1] [3] [6]. A prototype hybrid PV-biobattery Virtual Power Plant (VPP), consisting entirely of nanocomputers, microcontrollers, and high-frequency IoT sensors for microclimate, soil, and electrical parameters, was deployed and operated in Côte d’Ivoire to generate field validation datasets.
Reported P-MFC open-circuit voltages typically range from 0.2 V to 1.0 V per cell, yielding power densities between 1 and 100 mW⋅m−2 of anode area, subject to strong dependencies on electrode geometry, soil moisture, and plant physiology [1] [2] [7]. Existing modeling approaches generally pair electrochemical kinetics (Nernst equilibrium potentials combined with Butler-Volmer or Tafel expressions [8]) with porous-medium ionic transport, most notably via Archie’s law [9] or Monod-type microbial kinetics [10]. However, fully coupled, spatially resolved, transient multiphysics models capable of simultaneously resolving Richards-type unsaturated moisture transport [11], thermal transfer, ionic conduction, microbial substrate kinetics, and biofilm growth within a unified soil-electrode continuum remain underrepresented in the literature. Most existing works rely on lumped equivalent-circuit approximations or decoupled single-physics treatments.
This study investigates whether a physically derived, spatially resolved multiphysics framework can mechanistically capture and explain four distinct system behaviors:
1) the series/parallel electrical scaling under cell stacking;
2) the thermodynamic open-circuit potential of the selected bio-electrodes;
3) the dynamic operating point under real-world soil and biological loading;
4) the qualitative diurnal decoupling between a PV branch and a biopile branch operating in parallel.
This work evaluates the representative limits of the physics-based model across each regime using dedicated experimental validation datasets. The novel contributions of this work are threefold:
Analytical Derivation, a closed-form series-circuit solution for the coupled electrochemical sub-problem is derived from Kirchhoff’s laws and analytically proven to collapse into a direct dimensioning equation under explicit homogeneity and saturation bounds, providing a structural verification of the model topology.
High-Frequency Field Validation, the spatially resolved multiphysics model is confronted with high-frequency field current/voltage time-series (N > 6000 logs), moving beyond qualitative visual trend comparisons to quantify operational boundaries.
Mechanistic Decoupling Analysis, the distinct diurnal responses of the PV and bio-branches are reproduced purely through the differential structure of the coupled equations, providing a physical foundation for hybrid VPP control without empirical curve-fitting.
2. Equipment and Methods
2.1. Equipment
2.1.1. Modeling and Simulation Equipment
The modeling and simulation of the biopile and the photovoltaic system were performed using MATLAB software, which allows users to manipulate matrices and display graphs and data, among other things. MATLAB also enables users to implement algorithms and design user interfaces. It can communicate with various programming languages, including C, C++, Java, Python, and others.
2.1.2. Hardware Components of the Prototype Residential Virtual Power
Plant
Figure 1. Components of a residential VPP combining PV and a biofuel cell.
To compare the biosensor model’s simulation data with field data, a prototype virtual residential power plant based on the Internet of Things (IoT) was designed. It combines two energy sources: a solar panel and a plant-based biopile. This IoT device incorporates hardware components, including IoT sensors, microcontrollers, and a nano computer (Figure 1).
2.1.3. Software Used in the Residential Virtual Power Plant
On the software side, the VPP uses the Mosquitto MQTT broker to collect data acquired from each sensor via a microcontroller that contains an MQTT client. It uses Node-RED as an open-source visual programming tool based on data flows collected from the sensors, and PostgreSQL/PostGIS as a geographic database into which the sensor data is inserted after processing.
The NASA POWER web application was used to extract solar radiation data from the study area in order to compare it with data from the pyranometer installed on-site.
2.2. Methods
2.2.1. Physical and Mathematical Model of Energy Sources:
Biopile/Photovoltaics
1) Physical and mathematical model of the biopile
The biopile is modeled as a one-dimensional (1D) soil continuum of depth
, spatially discretized across
nodes spanning from the anode (
) to the cathode (
), with a nominal electrode cross-sectional area
(evaluated below as a calibrated parameter). Six spatially distributed state fields are resolved along this continuum, temperature
, volumetric water content
, dissolved oxygen concentration
, soluble substrate concentration
, proton concentration
, and planktonic microbial biomass
, coupled with three zero-dimensional (0D) lumped interfacial boundary states: anodic biofilm areal density
, residual magnesium anode quantity
, and passivation film thickness
. The complete dynamic system state vector
thus comprises
state variables.
a) Electrochemical sub-model
Anode and cathode equilibrium potentials are governed by the Nernst relation [8]:
(1)
(2)
where
and
represent the standard reduction potentials for the sacrificial magnesium/carbon electrode couple. In the low-overpotential regime, charge-transfer resistances (
) are formulated via the linearized Butler-Volmer equation [8]:
(3)
where the exchange current density
follows an Arrhenius temperature dependence with an activation energy of
. To account for bio-catalytic enhancement, the anodic exchange current density
is scaled by up to a factor of 5 as a function of the anodic biofilm density
, saturated via a half-saturation constant
. Ohmic transport within the unsaturated soil pore water is governed by Archie’s empirical law [9], expressed in conductivity form and modified by an Arrhenius thermal correction factor (
):
(4)
where
is the soil porosity,
is the volumetric water content,
is tortuosity, and
are empirical cementation and saturation exponents, respectively.
Given the strictly one-dimensional geometry and the lack of lateral current pathways, the quasi-static electrical sub-problem yields an exact closed-form solution derived from Kirchhoff’s Voltage Law (KVL) around the anode-electrolyte-cathode-Rext loop:
(5)
This closed-form formulation replaces prior iterative finite-volume numerical schemes, providing an exact analytical expression that collapses strictly into the algebraic series/parallel topology equations under uniform, steady-state, and saturated conditions.
b) Transport, kinetics, and thermal sub-models
Substrate transport (S) is governed by Fickian diffusion with an Arrhenius-corrected diffusion coefficient (DS) and a Monod-type microbial consumption reaction [10]:
(6)
incorporating a first-order biomass decay rate (
). To account for plant activity, rhizodeposition is introduced as a distributed substrate source within the root zone, spatially weighted by a Gaussian root-density distribution
centered at depth
with spatial variance
:
(7)
Dissolved oxygen (
) diffuses with diffusivity
, subject to a volumetric sink driven by aerobic substrate oxidation (
), a cathodic electrochemical flux boundary condition
, and atmospheric exchange
(
) at the top cathode surface. Proton dynamics (
) with diffusivity
are driven at the electrode boundaries by Faradaic anodic and cathodic current fluxes,
and
, respectively.
Soil moisture transport (
) is modeled via a diffusive Richards-type equation [11], subjected to a sinusoidal diurnal surface-humidity flux boundary condition (period 86,400 s) and a uniform root-extraction sink term
. At the anodic interface (x = 0), the growth of the attached biofilm areal density (
) follows a surface-restricted, Arrhenius-corrected Monod kinetic expression:
(8)
Thermal transport follows a 1D heat conduction equation parameterized by
thermal conductivity ktherm and volumetric heat capacity
, driven by internal Joule heating
, exothermic bio-oxidation heat generation
, and a surface boundary condition combining net absorbed solar radiation
and convective air cooling
. Finally, the depletion of the sacrificial magnesium anode follows Faraday’s law of electrolysis:
(9)
while the passivation layer dynamics are governed by a stoichiometric precipitation-dissolution mass balance.
c) Reduction to a lumped algebraic cell model
To verify the numerical implementation of the spatially resolved formulation, an asymptotic reduction to its lumped algebraic limit was performed. Under the four uniform and saturating hypotheses of a spatially uniform thermal field (
), constant soil moisture at reference saturation (
), non-limiting substrate supply (
), and a fully mature anodic biofilm (
), the distributed continuum field equations collapse to a single lumped operating point. Under these bounds, the spatial integral for internal resistance simplifies to scalar expressions for open-circuit voltage (
) and total internal resistance (
), governed directly by Equations (3)-(4). Direct numerical evaluation of the full 1D spatially resolved code under these exact hypotheses with a default load resistance
yields
,
,
, and
. These numerical outputs match the analytical scalar solution identically, confirming the structural consistency of the spatial discretization and ODE solver implementation by demonstrating that the continuum model preserves exact mathematical convergence to its lumped theoretical limit prior to experimental validation.
d) Parameter classification
Model parameters are systematically classified below according to their functional role and derivation origin, categorized as directly measured, literature-derived, or calibrated/adjusted values (Table 1). To maintain full methodological transparency, parameters designated as calibrated or adjusted are explicitly identified as fitted metrics and are not presented as independently measured physical properties.
Table 1. Classification of biopile parameters.
Parameter |
Symbol |
Value |
Unit |
Category |
Origin |
Column depth |
|
0.05 |
m |
Geometric
input |
Estimated
(device geometry) |
Electrode area |
|
1 × 10−3 (default);
0.18 - 0.2
(calibrated) |
m2 |
Geometric
input |
[Data Required]—not measured;
calibrated value is a curve fit, not an
independent measurement |
Anode standard potential |
|
−2.356 |
V |
Electrochemical |
Literature (Mg standard reduction
potential) |
Cathode standard potential |
|
0.401 |
V |
Electrochemical |
Literature (O2 reduction, alkaline) |
Anodic exchange
current density (ref.) |
|
0.05 |
A·m−2 |
Electrochemical |
Estimated (order-of-magnitude,
not fitted to this device) |
Cathodic exchange
current density (ref.) |
|
0.02 |
A·m−2 |
Electrochemical |
Estimated |
Pore-water conductivity |
|
0.5 |
S·m−1 |
Transport |
Estimated (typical moist-soil pore water) |
Porosity |
|
0.45 |
– |
Transport |
Estimated (typical loam) |
Archie
cementation/saturation exponents |
|
2.0, 2.0 |
– |
Transport |
Literature (Archie’s law, typical soil range) |
Max. specific growth
rate (ref.) |
|
1.5 × 10−5 |
s−1 |
Kinetic |
Estimated |
Monod half-saturation
constant |
|
0.5 |
mol·m−3
(model units) |
Kinetic |
Estimated |
External load resistance
(nominal script default) |
|
50 (default);
5, 7.66, 41.4 |
Ω |
Circuit input |
R_ext scenario values back-calculated
from measured V/I, not independently set |
Reference solar irradiance
amplitude (synthetic forcing) |
|
800/500 |
W·m−2 |
Input (climate forcing) |
Estimated (typical tropical clear-sky order
of magnitude); not fitted to the 15/10/2023 satellite record |
2) Photovoltaic branch model
The photovoltaic (PV) branch evaluated in the diurnal response analysis is modeled using a single-diode open-circuit voltage relationship, which captures the logarithmic dependence of the potential on solar irradiance (G) [12]:
(10)
where
denotes the number of series-connected solar cells,
is the diode ideality factor,
is the thermal voltage,
W⋅m−2 is the standard test condition (STC) irradiance, and
represents the nominal open-circuit voltage at STC. A hard lower bound (
) models the operational threshold below which the maximum power point tracking (MPPT) power electronics disengage.
This instantaneous algebraic formulation contrasts the fundamental time scales with the branch of the biosystem related to the physical and mathematical model, whose dynamic response is governed by biological and transport time constants spanning several days. This structural distinction provides the physical mechanism underlying the asymmetry of the diurnal response analyzed in the experimental analysis and dynamic simulation of a photovoltaic system.
2.2.2. Development of a Residential Virtual Power Plant
Combining Photovoltaic Technology and a Plant-Based Biopile
1) Hardware and software design of the residential Virtual Power Plant
The residential virtual power plant, built on an IoT architecture (Figure 2), consists of a photovoltaic component, a plant-based biopile component, and a battery.
Figure 2. Block diagram of the IoT system for remote monitoring of hybridization (PV, Biopile).
2) Design of the data Storage and processing system for the residential Virtual Power Plant
Figure 3. Data Storage and processing processes for the IoT remote monitoring device.
Data management is crucial in the residential virtual power plant with remote monitoring, as it enables the collection, storage, and processing of information from sensors. The process, from initial data acquisition to storage in a database, uses the MQTT protocol (Figure 3).
This general method shows that the sensors collect data on the renewable energy source, including temperature and energy production. This data is then structured into MQTT messages, each with a specific topic corresponding to the type of information to be transmitted (Table 2).
Table 2. Description of the publishing/subscription system.
Experiment |
Topic |
Description of the data collected |
Observation |
Energy
production |
INA1/production |
- Pot 1: Biopile Indicator or Control
without plant |
Current/Voltage
Measurement |
INA2/production |
- Biopile (single or series/parallel
combination) |
INA3/production |
- Photovoltaic (PV) panel |
INA4/production |
- Battery |
METEO1 |
- Pyranometer - Index pyranometer - Air humidity - Air temperature |
Weather
Parameter
Measurement |
METEO2 |
- Index temperature - Rain sensor - Light sensor |
SOL1 |
- Soil moisture sensor |
Soil Parameter
Measurement
for Plants |
SOL2 |
- NPK sensor: Potassium (N) - Nitrogen (P) - Phosphorus (K) |
Figure 4. Storage and processing of hybridization data (PV, Plant Biopile).
Having designed plant-based biofuel cells made of spirally arranged ginger with spikes and electrodes (magnesium and carbon), the next step is to produce the IoT remote monitoring system for the micro-generator. This system combines a photovoltaic panel and a plant-based microbial fuel cell to charge a battery and power electrical loads. During the production phase, all sensors send their messages to the MQTT server (Mosquitto), also known as the MQTT broker. Subscribers (clients) receive the MQTT messages relevant to their topics of interest and extract the useful data. This extracted data is then stored in a database, with tables structured by topic to facilitate information management. Node-RED acts as a software layer that simplifies the development, management, and automation of data flows and interactions between IoT devices in a remote monitoring system. It provides an efficient way to connect hardware and software components and offers a visual interface for creating and managing these flows (Figure 4).
3) Sizing of Biopiles
Figure 5. General method for sizing a plant-based biobattery microgenerator.
Solar energy and plant-based (phyto-microbial) fuel cells are renewable energy sources that are intermittent and depend on weather conditions, soil conditions, and plant cultivation conditions. Optimization ensures a continuous energy supply by adapting production to demand. Thus, this part of the research aims to contribute to scientific development by providing concrete applications for the design of plant-based biobatteries, taking into account both variations in solar energy and the changing conditions under which these biobatteries may be used.
By creating a mixed series/parallel configuration of a set of plant-based biosolar cells, it is possible to obtain a significant energy source capable of powering solar-powered electrical loads such as light bulbs, televisions, etc.
The method for optimizing plant-based biosolar cells is illustrated in Figure 5 above.
3. Results
3.1. Analysis of the Digital Implementation of the Physical and
Mathematical Model of the Biopile
All numerical results presented and developed in the mathematical modeling of the biobattery were generated using the custom multiphysics framework. The spatial discretization of the six dynamic fields distributed across the 21-node domain uses a second-order centered finite-difference scheme. The boundary conditions described in the biopile model are implemented using one-sided finite-difference closures, without ghost nodes, to ensure flux conservation at reactive spatial boundaries.
The resulting stiff system of 129 time-dependent ordinary differential equations (ODEs) was solved using the Radau implicit Runge-Kutta scheme (SciPy implementation) with relative and absolute integration tolerances set to 10−9 and 10−6, respectively. A parallel reference code was implemented using MATLAB’s ode15s stiff solver with variable step size and order for structural cross-validation. Both numerical implementations share the same basic equations, spatial discretizations, and parameter sets. Although the stability of the solution was maintained for all simulations, a formal study of mesh and tolerance independence remains a goal for future, more in-depth computational work.
The quasi-static electrical subproblem is solved algebraically at each step of the time evolution, directly from the instantaneous values of the spatial state. This direct algebraic evaluation eliminates numerical iteration errors and the risk of non-convergence. This formulation is based on the physical assumption that the electrical relaxation kinetics are decoupled from and significantly faster than the thermal, mass transport, and biological timescales—an assumption consistent with the estimated RC time constants of the system (Rint ≈ 1.5 kΩ) for transport dynamics spanning several hours and biological dynamics spanning several days.
3.2. Empirical Basis for Species and Electrode Selection
3.2.1. Empirical Basis for Selecting Plant Species for a Biopile
The classic microbial fuel cell transforms the chemical energy of organic matter into electricity by coupling with a living plant. This plant provides organic matter to the soil via photosynthesis and excretes substances through its roots, which soil microorganisms consume, thus releasing electrons collected by conductive electrodes [13].
The plant species were selected during a preliminary empirical screening campaign conducted on the physical prototype prior to the formulation of the model. Four candidate plant species were evaluated in parallel under identical instrumentation and monitoring conditions, yielding open-circuit potentials of 1.391 V for Solanum lycopersicum (tomato) and 1.643 V for Costus spicatus (spiral ginger), 1.352 V for Capsicum annuum (chili pepper), and 1.501 V for aloe vera (Figure 6). Costus spicatus exhibited the highest open-circuit response and was selected as the botanical matrix for all subsequent characterizations. Consequently, the rhizodeposition parameters (
,
,
) defined in the model are specific to the species Costus spicatus and do not represent generalized P-MFC values.
Figure 6. Measuring the energy produced by plants.
3.2.2. Empirical Bases for Biobile Electrodes
Anode and cathode materials are crucial for oxidation and reduction reactions. Anodes must maximize cell adhesion and electron transfer while being conductive, biocompatible, and chemically stable. Metallic materials used in anodes must be corrosion-resistant and non-toxic.
A subsequent screening evaluated four electrode pairs within Costus spicatus. Costus spicatus matrix: a magnesium anode paired with a graphitic carbon cathode (1.0 - 2.5 V), a magnesium anode paired with a charcoal cathode (1.0 - 1.5 V), a copper anode paired with a charcoal cathode (0.5 - 1.3 V), and an iron anode paired with a charcoal cathode (0.5 - 1.0 V) (Table 3). These voltage ranges are consistent with the values reported for comparable galvanic metal-carbon sacrificial pairs in ground energy recovery systems [1] [7]. The magnesium/graphitic carbon pair was selected for the final architecture to maximize the output potential, thereby forming the basis for the standard thermodynamic reduction potentials (
,
) presented in the biobattery model. The range of 1.0 - 2.5 V established during this preliminary screening serves as the empirical open-circuit reference evaluated in this study.
Table 3. Electrode configuration of a biobattery.
Selected Plant |
Performance Rating |
Electrodes |
Energy Efficiency |
Availability of electrodes |
Costus spicatus (spiked spiral ginger) |
Gold |
Anode: Magnesium Cathode: Carbon (Graphite) |
1 to 2.5 V |
Rare and must be ordered
from an outside supplier |
Silver |
Anode: Magnesium Cathode: Charcoal |
1 to 1.5 V |
Charcoal is available,
but magnesium is rare |
Bronze |
Anode: Copper (plumbing) Cathode: Charcoal |
0.5 to 1.3 V |
Available |
Platinum |
Anode: Iron Cathode: Charcoal |
0.5 to 1 V |
Widely available |
The measured voltage of a pot without a plant, using a substrate similar to that of the plant pots, was 0.916 V, while a pot containing sterilized sandy soil showed 0.840 V (Figure 7). These results were obtained under identical experimental conditions, including the same mass, water volume, electrode size, anode-cathode distance, and sunlight and temperature conditions.
Figure 7. Measurement of the energy produced by control pots without plants using the same Mg/graphite electrodes.
3.2.3. Residential Virtual Power Plant Prototyping
Figure 8 shows the residential virtual power plant (VPP) used to collect data via sensors and store it in a nanocomputer. The VPP consists of a photovoltaic solar panel, a biofuel cell array connected in parallel/series, a solar charge controller, a battery, and sensors concealed on-site in strategic locations, including planters, junction boxes, etc. This device is built on the architecture of the Internet of Things (IoT).
Figure 8. Presentation of the residential virtual power plant.
3.3. Analysis of Experimental Data Sources and Evaluation of
Model Performance
3.3.1. Analysis of Experimental Data Sources
The model’s performance was evaluated using two complementary empirical datasets collected from the physical prototype during a field campaign conducted on October 15, 2023, in Côte d’Ivoire.
The first source includes the architectural specifications and operating parameters defined during system assembly, prior to model development. These parameters include: 1) the nominal electrical specifications of the individual cells (0.5 - 1.5 V, operating range of 0.2 A), 2) the sizing diagrams for the stack topology for four series/parallel configurations (Figures 9-12), 3) open-circuit test parameters for the bioelectrodes (1.0 - 2.5 V for the selected magnesium/graphitic carbon pair), and 4) dashboard time series covering the photovoltaic (PV) voltage, biosensor voltage, battery status, ambient temperature, relative humidity, global solar irradiance, brightness, and precipitation.
In these log files, the PV branch shows nominal diurnal voltage fluctuations between 15.5 V and 20.0 V before experiencing a sharp drop at dusk starting at 18:57:50 and reaching 0.00 V at 19:25:16, whereas the potential of the biopile branch remains nearly stable at 2.50 V with a nominal current of 0.2 A throughout the recording, serving as the experimental open-circuit voltage
reference evaluated in the biopile model.
The second empirical source consists of raw, high-frequency digital IoT sensor logs exported directly from the data acquisition hardware (INA219 current/voltage/power monitors), synchronized with local soil/microclimate instrumentation and hourly satellite-derived solar irradiance records (NASA-POWER, ALLSKY_SFC_SW_DWN). To maintain complete methodological transparency, the spatial and parametric alignment between the exported raw spreadsheet channels and the rendered dashboard charts is treated as logically inferred.
Figure 9. Design of plant-based biopile: 4S3P.
Figure 10. Design of plant-based biopile: 8S3P.
Figure 11. Design of plant-based biobatteries: 6S2P.
Figure 12. Sizing of plant-based biopiles: 12 S2P.
3.3.2. Model Performance Evaluation
The model evaluation is organized according to three distinct levels of evidential rigor in order to explicitly distinguish between mathematical consistency, parametric calibration, and structural validation.
First, the reduction test performed within the model verifies internal mathematical consistency by comparing the 1D spatial implementation to its global algebraic limit, without relying on empirical data.
Second, the unmeasured geometric and circuit parameters (
and
) are calibrated directly from specific empirical operating points; the agreement observed in these cases represents a parametric fit rather than an independent validation and is referred to as such throughout the document.
Third, the circuit topology and the model’s emerging diurnal dynamics are evaluated against the system’s nominal specifications and raw telemetry records, without tuning the parameters to these specific target outputs, providing a structural evaluation limited to topological scaling and transient behavior. The crucial point is that the experimental data used for parametric calibration are not subsequently reused to claim independent validation of the same output.
A point-by-point statistical comparison of the time series (e.g., RMSE, Pearson correlation coefficient r) between the transient simulation and the raw sensor recordings was deliberately avoided due to the stationarity of the signal. The voltage/current recordings, whether single-cell or multi-cell, remain stationary at the sensor resolution throughout the 5.79-hour monitoring window, which makes the correlation coefficient statistically degenerate relative to a reference signal with near-zero variance. Consequently, the model’s performance is evaluated using the absolute and relative biases between the time-averaged empirical operating points and the model’s predictions, thereby providing a physically meaningful measure of error for stationary field data.
3.3.3. Loaded Operating Point: Single-Cell Branch
The raw data from the INA219 digital sensor for the single-cell monitoring channel (N = 6418 samples from the 5.79-hour field session, from 15:08:29 to 20:55:00) shows steady-state behavior, with an average cell voltage of V = 1.5025 ± 0.0076 V (range: 1.48 - 1.53 V) and an average output current of I = 0.1961 ± 0.0193 A, resulting in an implied external load resistance of
. This empirical voltage (1.5025 V) corroborates the nominal upper limit of the per-cell potential (1.5 V) used in this study, based on an independent, directly recorded dataset rather than a value transcribed from specification tables.
Figure 13 illustrates this empirical trajectory as well as the model’s response. The measured voltage and current profiles exhibit quasi-stationarity over the entire 5.79-hour monitoring period, indicating that the branch operated under a fixed external load, without dynamic disturbances. Physically, this invariant response is consistent with a static resistive load connected to an electrochemical matrix whose biological and transport dynamics evolve on characteristic time scales of the order of days rather than hours. A representative model output (
, corresponding to the upper limit of the parameter scan range) is plotted to visually illustrate the basic electrical behavior. Crucially, this graphical comparison highlights the two-order-of-magnitude discrepancy between the predictions based on default parameters (
) and the field measurements, thereby establishing the physical justification for the calibration sweep of the effective area.
![]()
![]()
Figure 13. Single-cell branch: empirical trajectory and model response.
The raw INA219 signal (blue, full 5.79-hour session) compared to the model prediction for
. and
, (red dashed line, an illustrative data point from the
sweep in the table below). The actual signal is virtually flat, given the instrument’s resolution; the model curve corresponds to one of the swept
values, not to a fitted trajectory. For a default electrode area
(
) under a fixed load
., the model predicts a quasi-steady-state output of I = 2.07 mA and V = 0.0158 V, underestimating the empirical operating point by two orders of magnitude due to an uncalibrated scaling of the internal resistance. To evaluate the model’s sensitivity and identify the effective active area required to reproduce field performance, a parametric sweep was performed while keeping all other kinetic and transport parameters constant (Table 4).
Table 4. Parametric scanning based on kinetic and transport parameters.
A_e (m2) |
I model (mA) |
V model (V) |
Bias (measured) ΔV (V) |
Ratio
I_model/I_measured |
1 × 10−3
(script default) |
2.07 |
0.016 |
−1.487 |
0.0105 |
1 × 10−2 |
19.96 |
0.153 |
−1.350 |
0.102 |
3 × 10−2 |
54.77 |
0.420 |
−1.083 |
0.279 |
5 × 10−2 |
83.97 |
0.643 |
−0.859 |
0.428 |
1 × 10−1 |
139.74 |
1.071 |
−0.432 |
0.713 |
2 × 10−1 |
209.09 |
1.602 |
+0.099 |
1.066 |
5 × 10−1 |
297.77 |
2.281 |
+0.779 |
1.518 |
The voltage bias (
) transitions through zero between
(
) and
(
), indicating that an effective active area of
(corresponding to an approximate 40 cm × 45 cm buried electrode surface) reconciles the model prediction with the empirical operating point. This result is presented explicitly as a parametric calibration rather than an independent structural validation, as
functions as a lumped scaling parameter that absorbs combined uncertainties from exchange current density (
) and soil electrolyte conductivity (
). Across the 500-fold variation in
(
to 0.50 m2), the predicted current scales nonlinearly by a factor of ~144 (from 2.07 mA to 297.77 mA), adhering directly to the theoretical inverse ohmic scaling law
and confirming that the underlying electrical formulation operates as derived regardless of the specific calibrated
value.
3.3.4. Circuit Topology: Series/Parallel Stacking Law
Applying Kirchhoff’s circuit laws to a stack comprising ns series-connected cells and np parallel branches yields the fundamental electrical scaling relation:
(11)
where
and
represent the nominal per-cell operating bounds established during preliminary characterization. Equation (11) calculates the nominal topological scaling according to Kirchhoff’s laws, compared to the design specifications of the 4S3P, 8S3P, 6S2P, and 12S2P configurations. This allows for verification of the consistency of the data in Table 5, which shows the measured and verified operating points for the 6S2P configuration against the designed model.
Table 5. Comparison of energy output between the model and the physical system.
Config. |
Cells |
V range
reported (V) |
ns× [0.5, 1.5] V, model (V) |
I reported (A) |
np× 0.2 A, model (A) |
Power
reported (W) |
Power
model (W) |
4S3P |
12 |
5 - 6 |
2.0 - 6.0 |
ND |
0.6 |
ND |
1.2 - 3.6 |
8S3P |
24 |
6 - 12 |
4.0 - 12.0 |
ND |
0.6 |
ND |
3.6 - 7.2 |
6S2P |
12 |
5 - 9 |
3.0 - 9.0 |
0.3 |
0.4 |
1.5 - 2.7 |
1.2 - 3.6 |
12S2P |
24 |
10 - 18 |
6.0 - 18.0 |
ND |
0.4 |
ND |
2.4 - 7.2 |
ND: Not Defined.
The 6S2P configuration data automatically recorded in the database during the production period made it possible to create Figure 14 below.
Figure 14. Experimental characterization of the electrical behavior of the biofuel cell_6S2P.
The maximum measured battery voltage corresponds precisely to ns × 1.5 V for all four configurations, indicating that the individual cells are operating near their upper thermodynamic limit in open circuit. This result is independently corroborated by the measurement of the field strength of a single cell (1.50 V). The measured battery currents correspond to the theoretical model np × 0.2 A for all four configurations. For the 6S2P array, a current of 0.3 A is mathematically consistent with the measured power and voltage limits: 8.5 V × 0.3 A = 2.55 W exactly reproduces the documented range of 2.0 to 3.6 W. This demonstrates that Kirchhoff’s scaling law accurately describes the array topology.
3.3.5. Experimental Analysis and Dynamic Simulation of a Photovoltaic
System and a Biopile
The Biopile_6S2P channel shows an average voltage of 7 ± 2 V and a constant current of 0.3000 A, both below the nominal limits of 9 V and 0.4 A. Channel 2 is a potentiometer combination, confirming the 6S2P topology. The single cell of the biopile and the PV branch were simulated, and their time stability is recorded in Table 6.
Table 6. Simulated temporal stability of a single biofuel cell and PV branch.
Quantity |
Start (15:08) |
End (20:52,
post-“sunset”) |
Relative change |
VPV (diode model) |
20.16 V |
0.00 V |
−100% |
Vbiopile (quasi-open-circuit) |
3.08522 V |
3.07104 V |
−0.460% |
Anode temperature T (0) |
298.15 K |
302.70 K |
+4.6 K |
Anodic biofilm density Xb |
1.00 × 10−4 kg·m−2 |
1.22 × 10−4 kg·m−2 |
+22% |
Remaining Mg (nMg) |
0.020000 mol |
0.020000 mol |
≈ 0% (near-open-circuit ⇒ negligible
Faradaic consumption) |
Figure 15. Dynamic divergence between the simulated photovoltaic voltage and the quasi-open-circuit potential of the biopile.
Figure 15 compares the simulated trajectory of a photovoltaic generator with the model’s near-open-circuit biofuel cell voltage. It highlights a qualitative difference in the simulated responses. The associated field recordings show a decrease in photovoltaic current and a near-steady biofuel cell voltage around 2.50 V, while the current tends towards zero. A stable voltage close to open circuit does not demonstrate available power under load after the drop in photovoltaic production; the data therefore only confirm voltage persistence, not a usable power supply complement.
Simulated day/night contrast: the PV model voltage drops to zero, while the modeled quasi-open-circuit biofuel cell voltage varies by approximately 0.46% over the interval (Table 6). This simulation does not demonstrate either the load output power or the available power of the biofuel cell.
The different response times result from the instantaneous relationship of the model’s PV diode and the slowness of biological and transport processes. This is a mechanistic deduction drawn from the model’s structure. Field telemetry shows voltage trends consistent with a qualitative day/night contrast, but because the biofuel cell channel is close to open circuit and its galvanic and microbial contributions are not dissociated, the measurements establish neither plant-induced causal complementarity nor usable power after sunset.
The mechanism governing this dynamic contrast is structural rather than artificial: VPV. It is expressed as an instantaneous algebraic function of the solar flux via the open-circuit relationship of an ideal diode, whereas Vbiopile depends on the local Nernstian equilibrium. In the model’s architecture, solar irradiance is coupled to the bioelectrochemical potential only indirectly, through a cascading sequence of biological and transport phenomena on a daily timescale (photosynthetic exudation → substrate diffusion → microbial kinetics and biofilm growth). These biochemical transitions are governed by reference kinetic rate constants (
et
) associated with characteristic relaxation times of the order of 105 to 106 s (approximately 1 to 10 days).
Over the simulated 5.73-hour interval, the model produces a rapid drop in photovoltaic voltage and a voltage variation in the biofuel cell of approximately 0.46% (Table 6). This behavior is simulated and does not constitute independent validation of the biofuel cell dynamics; therefore, the quasi-open voltage should not be interpreted as the power delivered.
The predicted quasi-open-circuit voltage of the biobattery (3.07 - 3.09 V) consistently exceeds the empirical quasi-open-circuit potential (2.50 V) by a difference of ΔV = 0.585 V. This potential difference stems from the ideal equilibrium formulation adopted in the model: the Nernst equation calculates an uncompensated thermodynamic ceiling that omits anode passivation phenomena at rest, parasitic losses at mixed potential, and non-ideal activity coefficients within the complex electrolytic matrix of the soil.
Figure 13 shows the empirical voltage and current recordings across the photovoltaic (PV) system on October 15, 2023, which serve as the direct experimental counterpart to the dynamic simulation trajectory. The telemetry data reveal two dominant physical regimes: sustained voltage oscillations in the afternoon, ranging from 15.5 to 20.0 V, followed by a rapid, nearly vertical drop in potential, which began at 6:57:50 p.m. and ceased at 7:25:16 p.m.
PV telemetry records voltage fluctuations in the afternoon, followed by a drop beginning around 18:57:50 (Figure 16). The sunset threshold used in the simulation was selected with reference to this observed transition. The comparison therefore documents the calibration basis of the forcing parameter; it is not an independent validation of the sunset timing or the PV model. Figure 16 establishes the empirical baseline against which the simulated temporal asymmetry (Figure 15) and the parameterized cutoff of solar forcing (sunset) are calibrated and structurally verified.
Figure 16. Experimental voltage and current measurements of the photovoltaic system.
The PV channel, which is the device’s PV monitoring channel, shows a voltage that fluctuates between approximately 15.5 and 20 V during the afternoon, then drops sharply starting at approximately 18:57:50, reaching 0 V around 19:25:16.
Figure 17. Experimental characterization of the electrical behavior of the biobattery.
Figure 17 shows the biofuel cell monitoring channel during the same session. The terminal voltage remains around 2.50 V while the current tends towards 0 A during the monitoring period (15:08:29 - 20:51:00), including the drop in photovoltaic voltage. This recording describes the persistence of a voltage close to that of an open circuit; it does not demonstrate the usable power of the biofuel cell after the drop in photovoltaic production.
Physically, this undisturbed profile confirms that the biobattery operated in a quasi-open-circuit condition (I ≈ 0 A), with the voltage measured across the terminals reflecting the actual, uncompensated bioelectrochemical equilibrium potential. Furthermore, it demonstrates that the characteristic relaxation times of the underlying biological and mass-transfer dynamics are sufficiently long such that a single diurnal solar transition does not induce any measurable disturbance in the cell’s potential. Finally, Figure 14 provides the main empirical evidence for day/night functional asymmetry, serving as a reference for establishing the structural—rather than parametric—validity of the differential-algebraic model formulation.
The biopile channel allows monitoring of the biopile in a physical device, showing a flat voltage of 2.5 V and a low current tending toward 0 A between 3:08 PM and 8:51 PM. This data constitutes major evidence of the qualitative day/night contrast, which is structurally reproduced by the simulation presented in this study.
3.3.6. Solar Forcing: Comparison between Model Assumptions, Actual
Irradiance Data, and Data on the Transition to Darkness
The synthetic solar forcing profile was compared to hourly satellite irradiance, on-site measured photovoltaic voltage, and optical ADC recordings from October 15, 2023. The pyranometer was shaded after approximately 3:30 PM and therefore did not serve as a quantitative reference for irradiance. Because the satellite data are hourly, they cannot resolve sub-hourly transient variations. The sunset threshold was determined from the on-site observed photovoltaic transition and is a calibrated input, not an independent validation.
Figure 18. Comparison of satellite-measured solar irradiance and analytical cosine forcing.
Figure 18 compares the hourly satellite irradiance to the analytical forcing profile. Their evening transition times differ by approximately one hour: the satellite data indicate zero at 6:00 PM, while the forcing profile uses 6:57:40 PM. The hourly satellite data cannot resolve this sub-hourly transition. Since the model threshold was selected from the photovoltaic transition observed on site, this comparison describes a calibrated input data and not an independent validation.
The satellite irradiance record indicates a diurnal peak of 649.6 W∙m−2 at 11:00:00, decaying through 267.9 W∙m−2, 143.4 W∙m−2, and 37.4 W∙m−2 at 15:00:00, 16:00:00, and 17:00:00, respectively, before dropping to zero at the 18:00:00 bin. Interpreted in isolation, this satellite time-series places the darkness transition approximately one hour prior to the model’s parameterized sunset (
, or 18:57:40).
The on-site photovoltaic voltage begins to decrease around 6:57:50 p.m., and the optical analog-to-digital converter’s reading reaches its lowest noise threshold between 7:23 p.m. and 7:30 p.m. These observations made it possible to determine the model’s sunset threshold at 6:57:40 p.m. Their agreement with this threshold constitutes calibration evidence, not an independent test. The spatial resolution of the zero irradiance interval at 6:00 p.m. in the satellite data is too coarse to assess the transition occurring in less than an hour.
On-site synchronized sensors provide the empirical basis used to select the sunset threshold. Since these same observations influenced the choice of boundary conditions, they cannot independently validate it. Hourly satellite data provide only a coarser contextual comparison.
4. Physical Interpretation and Discussion
The structural integrity of the circuit topology evaluation relies exclusively on Kirchhoff’s laws as applied to the terminal boundaries of each cell. Regardless of specific electrochemical or mass-transfer formulations, this result is the most robust in this study, as it does not depend on kinetic parameters that are subject to considerable physical uncertainty.
Conversely, the systematic overestimation of the open-circuit voltage—where the thermodynamic equilibrium potential (EOCV = 3.085 V) exceeds both the theoretical range of the electrodes (1.0 - 2.5 V) and the logarithm of the near-idle field (2.50 V) of a polarization of ΔV = 0.585 V, stems directly from the assumptions of ideal equilibrium inherent in the uncompensated Nernst equation. In physical magnesium/carbon (Mg/C) cells integrated into the soil, the actual terminal potentials are constrained by the self-passivation of the anode in open-circuit conditions (in accordance with the film thickness state variable
implemented in the state system, but not yet coupled to the anodic equilibrium potential
), by the parasitic kinetics at mixed potential resulting from competing reduction reactions at the anode interface, and by non-unit activity coefficients characteristic of porous soil solutions with high ionic concentration.
Furthermore, the significant underestimation of the charging current and the voltage across the terminals in the default geometric configuration stems from the inverse ohmic scaling relationship
. The nominal electrode area (
) is significantly smaller than the effective conductive area derived from field measurements. Although recalibration to
to 0.18 - 0.20 m2 corrects this discrepancy, this adjustment represents a single-parameter optimization. Consequently, isolating the individual physical contributions of the active surface area (
), the exchange current density (
), and the conductivity of soil interstitial water (
) to the internal resistance (
) remains a poorly formulated mathematical problem without independent empirical determination of at least one constitutive parameter.
The quasi-steady voltage/current recordings do not allow for the establishment of significant intra-session correlation statistics. The comparison at sunset is limited to a temporal comparison with the same in-situ photovoltaic recording used to select the threshold; it does not constitute independent validation.
The differential equations of the model imply different response times for the photovoltaic system and the biofuel cell. The voltage recorded in the biofuel cell remains quasi-steady during the decrease in photovoltaic production, but the channel is close to an open circuit, and the galvanic contribution of magnesium is not isolated. Therefore, the observation does not allow for the establishment of usable power or a causal effect related to the plant.
To place these operational parameters in the broader context of the energy literature, published studies on plant-based microbial fuel cells (P-MFCs) generally report open-circuit potentials ranging from 0.2 to 1.0 V per cell for sediment- and rhizosphere-based architectures using carbon or graphite cathodes, with power densities in the range of 1 - 100 mW∙m−2 of projected anode area. The measured operating point of the current system (V = 1.5025 V, I = 0.1961 A) places the terminal voltage near the upper limit of the ranges reported in the literature, while exhibiting an output current significantly higher than that of conventional bio-electrochemical systems.
This operational discrepancy is due to the incorporation of a magnesium sacrificial anode, which replaces the inert carbon/graphite materials typically used in research on plant-based microbial fuel cells (P-MFCs). The main anodic pathway introduces a dominant galvanic thermodynamic potential (due to corrosion), which overlaps with the bio-electrochemical electron transfer. Consequently, this system cannot be considered a direct equivalent of pure P-MFC systems that rely exclusively on microbial metabolic pathways. This distinction explains the high voltage and current values observed, though it does not constitute a validation of the model, as no data from the literature was used to adjust or calibrate the parameters in this study.
This operational divergence stems from the integration of a sacrificial magnesium anode, replacing the inert carbon/graphite materials typically used in plant-derived microbial fuel cell (P-MFC) research. The primary anodic pathway introduces a dominant galvanic thermodynamic potential (due to corrosion), which is superimposed on the bioelectrochemical electron transfer. Consequently, this system cannot be considered a direct equivalent of pure P-MFC systems relying exclusively on microbial metabolic pathways.
Published studies on plant-based microbial fuel cells (P-MFCs) generally report voltages and power densities for configurations different from that of this Mg/C assembly. The power measured at the terminals should not be compared here as a direct measure of microbial generation, because magnesium corrosion can generate a galvanic current. The available measurements do not quantify this contribution separately.
The sacrificial magnesium anode creates a galvanic pathway that is superimposed on any bioelectrochemical electron transfer. The measured signal must therefore be interpreted as the combined response of the entire Mg/C electrochemical system. Measurements taken in a pot without a plant, using the electrodes (Magnesium, Carbon), confirmed that microbial activity within the plant could produce a small amount of plant-derived electricity.
5. Implications and Perspectives
For the prototyping of residential virtual power plant (VPP), the structural alignments demonstrated during the scaling of circuit topology and diurnal dynamics validate the application of this modeling framework for sizing series/parallel cell configurations and for the qualitative assessment of the complementarity of PV and Biopile power. However, the current formulation does not yet allow for the quantitative prediction of power at the load operating point, as the electrode calibration discrepancy and the unresolved divergence in the Biopile_6S2P channel persist. Resolving these specific operational limitations requires targeted and scientifically justified follow-up steps.
Direct physical measurements of the anode and cathode active surfaces within a functional pot module are necessary to transform a single-parameter adjustment into a verifiable geometric parameter. Furthermore, electrochemical impedance spectroscopy (EIS), performed directly across the physical terminals of the stack (a standard characterization protocol in the literature on microbial and soil fuel cells), will allow us to separate the charge transfer resistance (
) from ohmic losses (
), thereby validating the Butler-Volmer and Archie submodels independently rather than collectively.
Conducting day/night transient simulations, driven directly by empirical irradiance time series with sub-minute resolution (
, will eliminate temporal discretization artifacts from hourly satellite records and allow for the evaluation of sunset dynamics at the native resolution of on-site optical sensors. Finally, extending the formulation of the anodic equilibrium potential (
) to include a passivation-dependent correction—derived from the existing film thickness state variable (
) will allow us to directly verify whether the film growth dynamics account for the systematic shift of +0.585 V in the open-circuit voltage.
6. Conclusions
This study defines the validity domain of a multiphysical model of a hybrid Mg/C biopile coupled to a photovoltaic panel. Theoretically, the spatial formulation reduces rigorously to its global algebraic limit (
). The topological submodel predicts a series/parallel stack with no deviation. At the single-cell level, the predicted charging voltage (1.50 V) matches the telemetric measurements to within three significant figures (1.5025 V, 0.1961 A). Furthermore, the differential-algebraic structure reproduces the day/night asymmetry, the collapse of the photovoltaic potential, and the stability of the biopile cell, the latter being considered an uncalibrated emergent property.
However, several mechanisms remain unproven hypotheses. The electrode surface calibration (Ae ≈ 0.18 - 0.20 m2) represents a single-parameter fit in the absence of direct geometric measurements. The overestimation of the open-circuit voltage (ΔV = +0.585 V, comparing EOCV = 3.085 V to the resting value of 2.50 V) stems from the ideal Nernst equation, which omits anode passivation, mixed potentials, and non-ideal ionic activity. Finally, the model’s sunset threshold was set at 18:57:40, based on the photovoltaic transition observed on-site, without quantitative irradiance measurements at the end of the day due to shading of the pyranometer after 3:30 PM. This agreement is calibration-dependent and does not allow for independent validation. The hourly satellite data, with irradiance canceled at 18:00:00, has too high a resolution to describe a one-hour transition, providing only a contextual comparison.
A complete quantitative validation requires direct physical measurements of the effective electrode area (Ae), decoupling of the charge transfer resistance
from the ohmic resistance
by electrochemical impedance spectroscopy (EIS [8] [14]), and a sub-hourly solar forcing simulation
. Furthermore, the power measured at the terminals cannot be equated to direct microbial production, as corrosion of the magnesium anode generates a galvanic current. Control tests without plants, carried out with the same Mg/C electrode pair, illustrate this.
The cell placed in sterilized sandy soil delivered 0.840 V, and that placed in unsterilized growing soil 0.916 V, compared to 1.5 V for the planted cell. The galvanic contribution would thus represent approximately 56% of the measured voltage. The difference of 0.076 V between sterilized soil and unsterilized soil (approximately 5%) suggests a contribution from soil microbial activity, while the difference of 0.584 V (approximately 39%) between the unsterilized control and the planted cell suggests a contribution associated with the plant and its rhizosphere. These results must be confirmed by repeated tests, carried out on the same substrate and under identical measurement conditions. In conclusion, this model offers a useful framework for the topological sizing of biofuel cells and for the qualitative analysis of their complementarity with photovoltaics in virtual power plants (VPP). On the other hand, in its current state, it does not make it possible to quantitatively predict the power delivered under load.
Author Contributions
Conceptualization, Kpegne K. H. & Aka B.; methodology, Koffi K. B.; software, Kone D.; validation, Mensah E. F., N’guessan K. R. & Lagobo Z. C.; formal analysis, Aka B.; investigation, Kpegne K. H.; resources, Kpegne K. H., Aka B., Mensah E. F. & Koffi K. B.; data curation, Kpegne K. H.; writing—original draft preparation, Kpegne K. H. & Kone D.; writing—review and editing, Kpegne K. H. & Kone D.; visualization, Mensah E. F. & N’guessan K. R.; supervision, Aka B. & Koffi K. B.; project administration, Kpegne K. H. All authors have read and agreed to the published version of the manuscript.