replace 3d array with 2d array in add_psfc_increment

[RussTreadon-NOAA]

Follow on to #57 with 2d array replacing 3d array in add_psfc_increment (file src/netcdf_io/calc_analysis.fd/inc2anl.f90)

Resolves #56

[RussTreadon-NOAA]

@CatherineThomas-NOAA , when you have time would you check to see if the changes in this PR to src/netcdf_io/calc_analysis.fd/inc2anl.f90 allow the C1152 case to successfully run to completion.

[CoryMartin-NOAA]

Changes look good, thanks @RussTreadon-NOAA don't worry I won't click merge until we have confirmation from @CatherineThomas-NOAA that everything works

[CatherineThomas-NOAA]

Thanks @RussTreadon-NOAA! The code looks good to me. I'll try it at C1152, though I will likely need to wait until WCOSS comes back later this week. We've been having issues on Hercules lately.

[RussTreadon-NOAA]

@CatherineThomas-NOAA . Thanks for the update. Let's hold this PR until we confirm that it works for C1152.

[CatherineThomas-NOAA]

I can confirm that this fix allows the C1152 case to run.

I'm thinking more about this code in general though. I believe that this code change is equivalent to what was already there, but I'm not sure if what was there is what we should be doing?

In the GSI minimization, we compute a increment for ps and then project it through the column:
delp(:,:,k) = pssm * (bk5(k)-bk5(k+1)) * r1000

In this utility, it is only pulling the increment from the lowest level delp to represent the surface. But the surface pressure should be impacted by increments to all levels above it. Do we instead need something like...

read 3D delp_inc
do loop k
ps_inc = ps_inc + delp_inc(k) / (bk5(k) - bk(k-1))
end do

[CoryMartin-NOAA]

@CatherineThomas-NOAA

The equation in GSI you provided:
delp(:,:,k) = pssm * (bk5(k)-bk5(k+1)) * r1000
let's simplify it (remove the r1000)
delp(:,:,k) = pssm * (bk5(k)-bk5(k+1))
If we assume pssm is the surface pressure increment, then solving for it gets us:
pssm = delp(k) / (bk5(k)-bk5(k+1))

I think what the GSI code is doing is saying "the surface pressure increment is applied as a scaling factor to all layers based on their thickness", so thus we should be able to get that increment back, right?

[CatherineThomas-NOAA]

@CoryMartin-NOAA: I agree, this process should be reversible.

Since a change like this would alter results and need more thorough investigation, how about we merge this PR as is so the C1152 g-w PR can advance, and then open a new issue to investigate the ps inc calculation? This PR would continue to replicate the behavior already used in our experiments.

[CoryMartin-NOAA]

@CatherineThomas-NOAA works for me