User: szilagyi Date: 2014/01/07 04:37 PM

Modified: /branches/tapir/src/ NullNews_GetNews.F90

Log: work on variable time-step compatibility and on diagnostic code

File Changes:

Directory: /branches/tapir/src/ ===============================

File [modified]: NullNews_GetNews.F90 Delta lines: +65 -15 =================================================================== --- branches/tapir/src/NullNews_GetNews.F90 2014-01-07 22:36:50 UTC (rev 23) +++ branches/tapir/src/NullNews_GetNews.F90 2014-01-07 22:37:24 UTC (rev 24) @@ -12,6 +12,7 @@ use NullNews_ScriUtil, only: NullNews_ResetInactiveRe use NullInterp use NullDecomp_SpinDecomp, only: SpinDecompFilter + use NullDecomp_IO

implicit none

@@ -24,20 +25,22 @@ CCTK_REAL, dimension (:,:,:), allocatable, save::& beta, omega, beta_u, u0

- CCTK_REAL, dimension(3) :: BUF logical, save :: FirstTime = .True. logical, save :: ThisProc = .false. integer, save :: loc_q, loc_p, myproc

CCTK_INT, dimension(1) :: Tarr - character(len=500) :: message

CCTK_INT :: ip + CCTK_REAL :: dt

DECLARE_CCTK_FUNCTIONS DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS

+ logical :: truncate + truncate = (IO_TruncateOutputFiles(cctkGH) .ne. 0) .and. FirstTime + ! call CCTK_INFO("Null News, Calculating News") if (FirstTime) then FirstTime = .false. @@ -82,9 +85,16 @@ endif myproc = CCTK_MyProc(cctkGH) endif + + ! assumes 'time' is time on the new level + time_of_news = 0.5*(null_time+null_time_p) + dt=null_time-null_time_p + if (use_linearized_omega .eq. 0) then - call NullNews_integ_omega (cctkGH, null_lsh, cctk_delta_time, tmp_rgfn, tmp_rgfs, Uo, Un, omegao, omegan,betao,betan) - call NullNews_integ_comega (cctkGH, null_lsh, cctk_delta_time, tmp_cgfn, tmp_cgfs, omegao, omegan, comegao, comegan) + call NullNews_integ_omega (cctkGH, null_lsh, dt, & + tmp_rgfn, tmp_rgfs, Uo, Un, omegao, omegan,betao,betan) + call NullNews_integ_comega (cctkGH, null_lsh, dt, & + tmp_cgfn, tmp_cgfs, omegao, omegan, comegao, comegan) else call NullNews_lin_omega (cctkGH, null_lsh, zeta, Jn, omegan, comegan) endif @@ -107,19 +117,56 @@ ! mid level values, from the characteristic evolution

J = 0.5d0 * (Jn + Jo) - J_u = (Jn - Jo) / cctk_delta_time + J_u = (Jn - Jo) / dt J_l = 0.5d0 * (Jn_l + Jo_l) - J_l_u = (Jn_l - Jo_l) / cctk_delta_time + J_l_u = (Jn_l - Jo_l) / dt beta = 0.5d0 * (betan + betao) cB = 0.5d0 * (cBn + cBo) U = 0.5d0 * (Un + Uo)

! For psi4 calculation - beta_u = (betan - betao) / cctk_delta_time - U_u = (Un - Uo) / cctk_delta_time + beta_u = (betan - betao) / dt + + U_u = (Un - Uo) / dt U_l = 0.5d0 * (Un_l + Uo_l) U_l_l = 0.5d0 * (Un_l_l + Uo_l_l)

+ + if (debug_output == 1) then + + call NullDecomp_WriteCoefsForComplexFunc('J_scri', cctkGH, truncate, & + null_lsh, 2, zeta, J(:,:,1), J(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForComplexFunc('Ju_scri', cctkGH, truncate, & + null_lsh, 2, zeta, J_u(:,:,1), J_u(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForComplexFunc('Jl_scri', cctkGH, truncate, & + null_lsh, 2, zeta, J_l(:,:,1), J_l(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForComplexFunc('Jlu_scri', cctkGH, truncate, & + null_lsh, 2, zeta, J_l_u(:,:,1), J_l_u(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForRealFunc('B_scri', cctkGH, truncate, & + null_lsh, zeta, beta(:,:,1), beta(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForRealFunc('Bu_scri', cctkGH, truncate, & + null_lsh, zeta, beta_u(:,:,1), beta_u(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForComplexFunc('U_scri', cctkGH, truncate, & + null_lsh, 1, zeta, U(:,:,1), U(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForComplexFunc('Uu_scri', cctkGH, truncate, & + null_lsh, 1, zeta, U_u(:,:,1), U_u(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForComplexFunc('Ul_scri', cctkGH, truncate, & + null_lsh, 1, zeta, U_l(:,:,1), U_l(:,:,2), time_of_news) + + call NullDecomp_WriteCoefsForComplexFunc('Ull_scri', cctkGH, truncate, & + null_lsh, 1, zeta, U_l_l(:,:,1), U_l_l(:,:,2), time_of_news) + + endif + + call NullInterp_d1 (eth_U(:,:,1), U(:,:,1), 1, 1) call NullInterp_d1 (eth_U(:,:,2), U(:,:,2), 1, 1)

@@ -154,9 +201,6 @@ call NullNews_ResetInactive(null_lsh, eth_J) call NullNews_ResetInactive(null_lsh, beth_J)

- ! assumes 'time' is time on the new level - time_of_news = cctk_time - .5d0 * cctk_delta_time - ! mid level values, from quantities computed above

sigmaJo = sigmaJn; sigmaKo = sigmaKn; sigmauo = sigmaun @@ -170,7 +214,7 @@

call NullNews_Psi4 (cctkGH, null_lsh, tmp_cgfn, tmp_cgfs, Jo, Jo_l,betao, cBo, Uo, Uo_l,& sigmaJn, sigmaKn, sigmaun,sigmarn, sigmarun, sigmauun,& - sigmaJo, sigmaKo, sigmauo,sigmaro, sigmaruo, sigmauuo,Psi4(:,:,1:2),cctk_delta_time) + sigmaJo, sigmaKo, sigmauo,sigmaro, sigmaruo, sigmauuo,Psi4(:,:,1:2),dt)

Psi4(:,:,1:2) = Psi4(:,:,1:2)/omega**3

@@ -185,7 +229,7 @@ endif

if (DEBUG_skip_BondiNews.eq.0) & - call NullNews_News2BondiNews(cctkGH, null_lsh, cctk_delta_time, qsize,& + call NullNews_News2BondiNews(cctkGH, null_lsh, dt, qsize,& Uo, Un, Uyo, Uyn, zEvolo, zEvoln,& zEvolh, News(:,:,1:2), NewsB(:,:,1:2), uBondio, uBondi, uBondin,& redshiftB, omega, J, eth_J, beth_J, eth_U, beth_U,& @@ -194,8 +238,15 @@ n_tmp3_cgfn, n_tmp3_cgfs, n_tmp1_rgfn, n_tmp1_rgfs,& n_tmp2_rgfn, n_tmp2_rgfs, starting, EV_mask)

- dMdOmega = NewsB(:,:,1:2) * conjg(NewsB(:,:,1:2)) * redshiftB * circle + if (debug_output == 1) then

+ call NullDecomp_WriteCoefsForRealFunc('uB_scri', cctkGH, truncate, & + null_lsh, zeta, uBondi(:,:,1), uBondi(:,:,2), time_of_news) + + end if + + dMdOmega = real(NewsB(:,:,1:2) * conjg(NewsB(:,:,1:2))) * redshiftB * circle + if (mask_NewsB .ne. 0) then do ip = 1, 2 NewsB(:,:,ip) = EQ_mask(:,:) * NewsB(:,:,ip) ! EV_mask @@ -226,7 +277,6 @@ omegao = omegan comegao = comegan

- end subroutine NullNews_GetNews