  The energy functional, the governing partial differential equation(s) (PDE),
and the boundary conditions need to be consistent with each other in a modeling
system. In electrolyte solution study, people usually use a free energy form of
an infinite domain system (with vanishing potential boundary condition) and the
derived PDE(s) for analysis and computing. However, in many real systems and/or
numerical computing, the objective domain is finite, and people still use the
similar energy form, PDE(s) but with different boundary conditions, which may
cause inconsistency. In this work, (1) we present a mean field free energy
functional for electrolyte solution within a finite domain with either physical
or numerically required artificial boundary. Apart from the conventional energy
components (electrostatic potential energy, ideal gas entropy term and chemical
potential term), new boundary interaction terms are added for both Neumann and
Dirichlet boundary conditions. These new terms count for physical interactions
with the boundary (for real boundary) or the environment influence on the
computational domain system (for non-physical but numerically designed
boundary). (2) The traditional physical-based Poisson-Boltzmann (PB) equation
and Poisson-Nernst-Planck (PNP) equations are proved to be consistent with the
new free energy form, and different boundary conditions can be applied. (3) In
particular, for inhomogeneous electrolyte with ionic concentration-dependent
dielectric permittivity, we derive the generalized Boltzmann distribution
(thereby the generalized PB equation) for equilibrium case, and the generalized
PNP equations for non-equilibrium case, under different boundary conditions.
Numerical tests are performed to demonstrate the different consequences
resulted from different energy forms and their derived PDE(s).
