periodic2_ewald_exp_erfc Function

public pure function periodic2_ewald_exp_erfc(kz, karg0, adz) result(term)

periodic2 Ewald の k≠0 項 exp(kz)erfc(karg0 + adz) を overflow させずに返す。 kz = ±|k|dz と adz = ±alphadz は同じ符号で、karg0 = |k|/(2 alpha) なので kz - (karg0 + adz)^2 = -(karg0^2 + adz^2) が成り立つ。 exp(kz) が binary64 で有限な kz では元の積をそのまま返す(従来の式とビット一致)。 exp(kz) が overflow する kz では y = karg0 + adz >= sqrt(2 kz) > 37 なので erfc は 0 に underflow し、 元の積は inf0 = NaN になる。この範囲だけ exp(-(karg0^2 + adz^2))erfc_scaled(y) で評価する。 真値はこのとき exp(-kz) 以下で、binary64 では 0 か非正規数になる。

Arguments

Type IntentOptional Attributes Name
real(kind=dp), intent(in) :: kz

指数の引数 ±|k|dz。

real(kind=dp), intent(in) :: karg0

指数の引数 ±|k|dz。 |k|/(2 alpha)。

real(kind=dp), intent(in) :: adz

指数の引数 ±|k|dz。 |k|/(2 alpha)。 kz と同符号の ±alpha*dz。

Return Value real(kind=dp)


Called by

proc~~periodic2_ewald_exp_erfc~~CalledByGraph proc~periodic2_ewald_exp_erfc periodic2_ewald_exp_erfc proc~add_periodic2_exact_ewald_correction_single_source add_periodic2_exact_ewald_correction_single_source proc~add_periodic2_exact_ewald_correction_single_source->proc~periodic2_ewald_exp_erfc proc~add_periodic2_exact_ewald_field_potential_single_source add_periodic2_exact_ewald_field_potential_single_source proc~add_periodic2_exact_ewald_field_potential_single_source->proc~periodic2_ewald_exp_erfc proc~add_periodic2_exact_ewald_kneq0_correction_single_source add_periodic2_exact_ewald_kneq0_correction_single_source proc~add_periodic2_exact_ewald_kneq0_correction_single_source->proc~periodic2_ewald_exp_erfc proc~add_periodic2_exact_ewald_potential_correction_single_source add_periodic2_exact_ewald_potential_correction_single_source proc~add_periodic2_exact_ewald_potential_correction_single_source->proc~periodic2_ewald_exp_erfc proc~add_periodic2_kneq0_potential_correction_single_source add_periodic2_kneq0_potential_correction_single_source proc~add_periodic2_kneq0_potential_correction_single_source->proc~periodic2_ewald_exp_erfc proc~precompute_periodic_root_operator precompute_periodic_root_operator proc~precompute_periodic_root_operator->proc~add_periodic2_exact_ewald_field_potential_single_source proc~core_complete_plan_impl core_complete_plan_impl proc~core_complete_plan_impl->proc~precompute_periodic_root_operator proc~core_build_panel_plan_impl core_build_panel_plan_impl proc~core_build_panel_plan_impl->proc~core_complete_plan_impl proc~core_build_plan_impl core_build_plan_impl proc~core_build_plan_impl->proc~core_complete_plan_impl

Source Code

  pure real(dp) function periodic2_ewald_exp_erfc(kz, karg0, adz) result(term)
    real(dp), intent(in) :: kz, karg0, adz
    real(dp), parameter :: max_finite_exp_arg = log(huge(1.0_dp))

    if (kz <= max_finite_exp_arg) then
      term = exp(kz)*erfc(karg0 + adz)
    else
      term = exp(-(karg0*karg0 + adz*adz))*erfc_scaled(karg0 + adz)
    end if
  end function periodic2_ewald_exp_erfc