Add ISOLATE_MASS_WEIGHT_PGF option to PGF calculation

[claireyung]

This commit adds the runtime parameter ISOLATE_MASS_WEIGHT_PGF with
default value False. If true, then the PGF calculation is modified
so that MASS_WEIGHT_IN_PRESSURE_GRADIENT cannot affect cells above and
below through the summation of intx_dpa or intx_dza (and y equivalents).
This summation can cause problems near the grounding line of ice shelves,
even with RESET_INTXPA_INTEGRAL (which was a workaround for this problem),
if there is no suitable grid cell in the entire column that is
non-vanished and non-tilted. This code recalculates intx_dpa without
MASS_WEIGHT_IN_PRESSURE_GRADIENT and uses that non-weighted version for
the calculation of intx_pa in subsequent cells,but still
retains intx_dpa potentially affected by MASS_WEIGHT_IN_PRESSURE_GRADIENT
in the PFu and PFv calculation to retain its desired effect near vanished,
sloped layers susceptible to grid-scale noise.

The result is that in the ISOMIP+ test case, combined with #810
velocities near the grounding line are reduced to 10^(-12)m/s.

Currently the code is very inefficient. As of now, if ISOLATE_MASS_WEIGHT_PGF
= True, and MWIPG = True, it calculates ALL the density integrals twice (one with
and one without MWIPG). This involves lots of eqn of state calls so slows the
model down a lot. This should be improved (though I don’t have a better
idea right now).

Default answers should not change since default is False, keeping original
intx_dpa and allowing MWIPG in one cell to affect cells above/below.

I would appreciate suggestions to optimise this code!!