  Density dependence of nuclear symmetry energy as well as its partial wave
decomposition is studied within the framework of lowest-order constrained
variational (LOCV) method using AV18 two-body interaction supplemented by UIX
three-body force. The main focus of the present work is to introduce a revised
version of three-body force which is based on an isospin-dependent
parametrization of coefficients in the UIX force, in order to overcome the
inability to produce correct saturation-point parameters} in the framework of
LOCV method. We find that employing the new model of {\ph three-body force} in
the LOCV formalism leads to successfully reproducing the semi-empirical
parameters of cold nuclear matter, including} $E_{sym}(\rho_0)$, $L$, and
$K_{sym}$. All our models of three-body force combined with AV18 two-body force
give maximum neutron star mass higher than $2\;M_\odot$. The fraction of
protons in the nucleon cores of neutron stars strongly depends on the
three-body force parametrization.
