Water Bottle Rocket Thrust - two calculation methods not matching

The first method is correct. In the second you have assumed that the pressure at the nozzle is still $P$ despite the water exiting with some velocity. i.e. You have neglected the dynamic pressure.

You need to use Bernoulli's principle $$ P + \frac{\rho v^2}{2} + \rho h g = {\rm constant}$$

Your first method assumes that the top surface of the water is hardly moving (because it's surface area is much bigger than the nozzle area). Applying the same idea to the second method, then we can calculate the constant both in the water and immediately below the nozzle as $$P = P_A + \frac{\rho v^2}{2},$$ where $P_A$ is atmospheric pressure and we neglect the small $\rho h g$ term which increases pressure due to to the column of liquid above the nozzle on the LHS. If we further assume that $P \gg P_A$ then $P = \rho v^2/2$ and the rate of change of momentum of liquid from the nozzle is $$ F = \rho A v^2 = 2PA$$

Often these rockets are launched off a vertical section of pipe or rod that extends up inside the bottle, acting as a piston until the rocket has moved far enough to clear the end of the pipe. During this phase, your second method is correct: the thrust is simply the nozzle area times the pressure.

Why should the thrust increase when the rocket separates from the piston? Let me try to provide an intuitive justification. I won't prove that the thrust doubles, just dispel the notion that it should remain unchanged.

Let's say that somehow the piston is able to extend itself by constantly adding little cylindrical plugs to itself. These plugs are initially arranged in a tall rack just beside the rocket; as the rocket moves upward, the piston keeps grabbing a plug from the rack, somehow transporting it through the wall of the rocket, and adding it on to the end. The rocket never leaves the piston, and it's clear that the thrust remains just $PA$.

But this is essentially what is actually happening as the rocket expels water, with one exception. Each little bit of water leaving the nozzle can be thought of as a "plug", and the force acting to separate it from the rest of the rocket is still $PA$. But, unlike the plugs waiting in the rack at rest relative to the earth, each one of the water plugs is moving upward a bit faster than the previous one was--just fast enough to match the velocity of the rocket. The momentum transferred to the rocket by these moving plugs constitutes an additional force relative to the case in which the plugs have zero velocity.

Of course, this analogy of plugs being added during flight does not correspond to the actual case of a rocket that loses mass over time. But that difference does not affect the instantaneous thrust.

The pressure drop from tank pressure to atmospheric pressure does not occur instantly at the nozzle but rather is spread out according to the area of the flow channel. This reduced pressure results in additional thrust that was not accounted for in your second solution. Here is one way to calculate this additional missing thrust:

Bernoulli's equation (your second to last equation) applies:

$$P+\frac12\rho v^2 = constant$$

we can combine this with your mass flow equation to get:

$$P+\frac12\rho \left(\frac{\dot m}{\rho A}\right)^2 = constant$$

Your original answer assumes negligible velocity/large area at the water surface¹:

$$P_{tank} + 0 = constant$$

This gives us our constant:

$$P+\frac12\rho \left(\frac{\dot m}{\rho A}\right)^2 = P_{tank}$$

and we know the pressure at the exit is atmospheric/ 0 gauge pressure:

$$0 + \frac12\rho \left(\frac{\dot m}{\rho A_{exit}}\right)^2 = P_{tank}$$

We can solve for $\dot m$:

$$\dot m = A_{exit}\sqrt{2 \rho P_{tank}}$$

Plugging back in:

$$P+P_{tank}\left(\frac{A_{exit}}{A}\right)^2 = P_{tank}$$

Solving for pressure:

$$P = P_{tank} \left(1-\left(\frac{A_{exit}}{A}\right)^2\right)$$

So if we wanted to calculate the additional thrust due to the pressure ebeing lower near the nozzle we'd need to integrate the pressure by the area:

$$F = F_{up} - F_{down} = \int_{A_{exit}}^\infty \left (P_{tank} - P_{tank} \left(1-\left(\frac{A_{exit}}{A}\right)^2\right) \right) \; dA $$

$$F = P_{tank} \int_{A_{exit}}^\infty \frac{{A_{exit}}^2}{A^2} \; dA $$

$$F = P_{tank} A_{exit} $$

So there's the missing $PA$ from your second solution.

^{1: You can make your equations more accurate (especially for thin bottle rockets) by using the actual cross sectional area at the water surface instead of infinity both here and in the limit of the integral.}