diff options
-rw-r--r-- | src/FuncAndJacobian.c | 18 | ||||
-rw-r--r-- | src/TwoPunctures.c | 31 |
2 files changed, 26 insertions, 23 deletions
diff --git a/src/FuncAndJacobian.c b/src/FuncAndJacobian.c index 6f1b993..b2f265e 100644 --- a/src/FuncAndJacobian.c +++ b/src/FuncAndJacobian.c @@ -326,6 +326,14 @@ F_of_v (CCTK_POINTER_TO_CONST cctkGH, } if (do_residuum_debug_output) { + indx = Index (0, i, j, 0, nvar, n1, n2, n3); + al = Pih * (2 * i + 1) / n1; + A = -cos (al); + be = Pih * (2 * j + 1) / n2; + B = -cos (be); + AB_To_XR (nvar, A, B, &X, &R, U); + C_To_c (nvar, X, R, &x, &r, U); + rx3_To_xyz (nvar, x, r, 0.0, &y, &z, U); r_plus = sqrt ((x - par_b) * (x - par_b) + y * y + z * z); r_minus = sqrt ((x + par_b) * (x + par_b) + y * y + z * z); psi = @@ -334,11 +342,11 @@ F_of_v (CCTK_POINTER_TO_CONST cctkGH, psi4 = psi2 * psi2; psi7 = psi * psi2 * psi4; fprintf(debugfile, "%.8g %.8g %.8g %.8g %.8g %.8g\n", x, y, A, B, - U.d11[0] + - U.d22[0] + - U.d33[0] - + 2.0 * Pi / psi2/psi * sources[Index(0,i,j,k,1,n1,n2,n3)], - U.d0[0] + U.d11[indx] + + U.d22[indx] + + U.d33[indx] + + 2.0 * Pi / psi2/psi * sources[indx], + U.d0[indx] ); } } diff --git a/src/TwoPunctures.c b/src/TwoPunctures.c index 5d1f792..c3bd889 100644 --- a/src/TwoPunctures.c +++ b/src/TwoPunctures.c @@ -59,32 +59,28 @@ void set_initial_guess(CCTK_POINTER_TO_CONST cctkGH, rx3_To_xyz (nvar, s_x[i3D], r, phi, &(s_y[i3D]), &(s_z[i3D]), U); } Set_Initial_Guess_for_u(cctkGH, n1*n2*n3, v.d0, s_x, s_y, s_z); - - if (do_initial_debug_output) - debug_file=fopen("initial.dat", "w"); for (i = 0; i < n1; i++) for (j = 0; j < n2; j++) - { for (k = 0; k < n3; k++) { - indx=Index(0,i,j,k,1,n1,n2,n3); - tmp_r= s_x[indx]*s_x[indx]+ - s_y[indx]*s_y[indx]+ - s_z[indx]*s_z[indx]; - v.d0[indx]/=-tmp_r/((sqrt(tmp_r)+10)*(sqrt(tmp_r)+10))+1; + indx = Index(0,i,j,k,1,n1,n2,n3); + v.d0[indx]/=(-cos(Pih * (2 * i + 1) / n1)-1.0); } - indx=Index(0,i,j,0,1,n1,n2,n3); - if (do_initial_debug_output) + Derivatives_AB3 (nvar, n1, n2, n3, v); + if (do_initial_debug_output) + { + debug_file=fopen("initial.dat", "w"); + for (i = 0; i < n1; i++) + for (j = 0; j < n2; j++) { - fprintf(debug_file, "%.8g %.8g %.8g %.8g %.8g\n", s_x[indx], s_y[indx], + fprintf(debug_file, "%.8g %.8g %.8g %.8g %.8g\n", + s_x[indx], s_y[indx], v.d0[indx], - v.d0[indx]*(-cos(Pih * (2 * i + 1) / n1)-1.0), - (-cos(Pih * (2 * i + 1) / n1)-1.0)); + (-cos(Pih * (2 * i + 1) / n1)-1.0), + v.d1[indx]); } - } - if (do_initial_debug_output) fclose(debug_file); - Derivatives_AB3 (nvar, n1, n2, n3, v); + } free(s_z); free(s_y); free(s_x); @@ -124,7 +120,6 @@ TwoPunctures (CCTK_ARGUMENTS) if (use_external_initial_guess) { set_initial_guess(cctkGH, v); - F_of_v (cctkGH, nvar, n1, n2, n3, v, F, u); } Newton (cctkGH, nvar, n1, n2, n3, v, Newton_tol, Newton_maxit); |