Make `normal_lcdf(x | 0, 1)` equivalent to `Phi()`

Since https://github.com/stan-dev/math/pull/1411 normal_lcdf is more robust in the tails than Phi or even Phi_approx. Update Phi() to be as robust as normal_lcdf and add log_Phi() as discussed in the linked PR above (should this be a separate issue?)

In fact, Phi_approx() is less robust than Phi(). Adding @PhilClemson @nhuurre @bbbales2

generated quantities {
  vector[N] b;
  vector[N] c;
  vector[N] d;
   for (i in 1:N) {
    b[i] = normal_lcdf(i | 0, 1);
    c[i] = log(Phi_approx(i));
    d[i] = log(Phi(i));
   }
}

> data.table(mean_b = test_b, mean_c = test_c, mean_d = test_d)
           mean_b       mean_c       mean_d
 1:  -1.72754e-01 -1.72771e-01 -1.72754e-01
 2:  -2.30129e-02 -2.30241e-02 -2.30129e-02
 3:  -1.35081e-03 -1.23271e-03 -1.35081e-03
 4:  -3.16717e-05 -1.83432e-05 -3.16717e-05
 5:  -2.86652e-07 -5.01624e-08 -2.86652e-07
 6:  -9.86588e-10 -1.65181e-11 -9.86588e-10
 7:  -1.27981e-12 -4.44089e-16 -1.27987e-12
 8:  -6.22096e-16  0.00000e+00 -6.66134e-16
 9:  -1.12859e-19  0.00000e+00  0.00000e+00
10:  -7.61985e-24  0.00000e+00  0.00000e+00
11:  -1.91066e-28  0.00000e+00  0.00000e+00
12:  -1.77648e-33  0.00000e+00  0.00000e+00
13:  -6.11716e-39  0.00000e+00  0.00000e+00
14:  -7.79354e-45  0.00000e+00  0.00000e+00
15:  -3.67097e-51  0.00000e+00  0.00000e+00
16:  -6.38875e-58  0.00000e+00  0.00000e+00
17:  -4.10600e-65  0.00000e+00  0.00000e+00
18:  -9.74095e-73  0.00000e+00  0.00000e+00
19:  -8.52722e-81  0.00000e+00  0.00000e+00
20:  -2.75362e-89  0.00000e+00  0.00000e+00
21:  -3.27928e-98  0.00000e+00  0.00000e+00
22: -1.43989e-107  0.00000e+00  0.00000e+00
23: -2.33064e-117  0.00000e+00  0.00000e+00
24: -1.39039e-127  0.00000e+00  0.00000e+00
25: -3.05670e-138  0.00000e+00  0.00000e+00
26: -2.47606e-149  0.00000e+00  0.00000e+00
27: -7.38948e-161  0.00000e+00  0.00000e+00
28: -8.12387e-173  0.00000e+00  0.00000e+00
29: -3.28979e-185  0.00000e+00  0.00000e+00
30: -4.90671e-198  0.00000e+00  0.00000e+00
31: -2.69525e-211  0.00000e+00  0.00000e+00
32: -5.45208e-225  0.00000e+00  0.00000e+00
33: -4.06119e-239  0.00000e+00  0.00000e+00
34: -1.11390e-253  0.00000e+00  0.00000e+00
35: -1.12491e-268  0.00000e+00  0.00000e+00
36: -4.18262e-284  0.00000e+00  0.00000e+00
37: -5.72557e-300  0.00000e+00  0.00000e+00
38: -2.88543e-316  0.00000e+00  0.00000e+00
39:   0.00000e+00  0.00000e+00  0.00000e+00
40:   0.00000e+00  0.00000e+00  0.00000e+00

Start by reading the linked PR 1411 and the existing implementations of normal_lcdf, Phi, and Phi_approx. Use the generated-quantities comparison in the issue to check tail behavior; done means Phi matches the robust normal_lcdf behavior and log_Phi is added as discussed.