fmm_m2l_rotation_work Subroutine

private subroutine fmm_m2l_rotation_work(c, pm, pl, vscales, m2l_ztranslate_coef, alpha, src_m, beta, dst_l, work)

Direct M2L translation by 4 rotations and 1 translation

Compute the following matrix-vector product: \f[ \mathrm{dst} = \beta \mathrm{dst} + \alpha L_M \mathrm{src}, \f] where \f$ \mathrm{dst} \f$ is a vector of coefficients of output spherical harmonics, \f$ \mathrm{src} \f$ is a vector of coefficients of input spherical harmonics and \f$ L_M \f$ is a matrix of a multipole-to-local translation.

Rotates around OZ and OY axes, translates over OZ and then rotates back around OY and OZ axes.

@param[in] c: Radius-vector from new to old centers of harmonics @param[in] src_r: Radius of old harmonics @param[in] dst_r: Radius of new harmonics @param[in] pm: Maximal degree of multipole spherical harmonics @param[in] pl: Maximal degree of local spherical harmonics @param[in] vscales: Normalization constants for Y_lm @param[in] m2l_ztranslate_coef: @param[in] alpha: Scalar multiplier for src_m @param[in] src_m: Expansion in old harmonics @param[in] beta: Scalar multiplier for dst_l @param[inout] dst_l: Expansion in new harmonics @param[out] work: Temporary workspace of a size 6pp+19*p+8 where p is a maximum of pm and pl

Arguments

Type IntentOptional Attributes Name
real(kind=rp), intent(in) :: c(3)
integer, intent(in) :: pm
integer, intent(in) :: pl
real(kind=rp), intent(in) :: vscales((pm+pl+1)**2)
real(kind=rp), intent(in) :: m2l_ztranslate_coef(pm+1,pl+1,pl+1)
real(kind=rp), intent(in) :: alpha
real(kind=rp), intent(in) :: src_m((pm+1)*(pm+1))
real(kind=rp), intent(in) :: beta
real(kind=rp), intent(inout) :: dst_l((pl+1)*(pl+1))
real(kind=rp), intent(out), target :: work(6*max(pm,pl)**2+19*max(pm,pl)+8)

Calls

proc~~fmm_m2l_rotation_work~~CallsGraph proc~fmm_m2l_rotation_work fmm_m2l_rotation_work proc~carttosph carttosph proc~fmm_m2l_rotation_work->proc~carttosph proc~fmm_m2l_ztranslate_work fmm_m2l_ztranslate_work proc~fmm_m2l_rotation_work->proc~fmm_m2l_ztranslate_work proc~trgev trgev proc~fmm_m2l_rotation_work->proc~trgev proc~fmm_sph_rotate_oxz_work fmm_sph_rotate_oxz_work proc~fmm_m2l_rotation_work->proc~fmm_sph_rotate_oxz_work proc~fmm_sph_rotate_oz_adj_work fmm_sph_rotate_oz_adj_work proc~fmm_m2l_rotation_work->proc~fmm_sph_rotate_oz_adj_work proc~fmm_sph_rotate_oz_work fmm_sph_rotate_oz_work proc~fmm_m2l_rotation_work->proc~fmm_sph_rotate_oz_work

Called by

proc~~fmm_m2l_rotation_work~~CalledByGraph proc~fmm_m2l_rotation_work fmm_m2l_rotation_work proc~fmm_m2l fmm_m2l proc~fmm_m2l->proc~fmm_m2l_rotation_work proc~tree_m2l tree_m2l proc~tree_m2l->proc~fmm_m2l proc~cart_propnear_at_ipart cart_propnear_at_ipart proc~cart_propnear_at_ipart->proc~fmm_m2l proc~fmm_solve_for_multipoles fmm_solve_for_multipoles proc~fmm_solve_for_multipoles->proc~tree_m2l proc~fmm_solve fmm_solve proc~fmm_solve->proc~tree_m2l proc~cart_prop_at_ipart cart_prop_at_ipart proc~cart_prop_at_ipart->proc~cart_propnear_at_ipart proc~prepare_fmm_ext_ipd prepare_fmm_ext_ipd proc~prepare_fmm_ext_ipd->proc~fmm_solve_for_multipoles proc~preapare_fmm_static preapare_fmm_static proc~preapare_fmm_static->proc~fmm_solve_for_multipoles proc~field_extd2d field_extD2D proc~field_extd2d->proc~prepare_fmm_ext_ipd proc~prepare_fmm_ipd prepare_fmm_ipd proc~prepare_fmm_ipd->proc~prepare_fmm_ext_ipd proc~elec_prop_m2m elec_prop_M2M proc~elec_prop_m2m->proc~preapare_fmm_static proc~elec_prop_m2d elec_prop_M2D proc~elec_prop_m2d->proc~preapare_fmm_static proc~tmatvec_otf TMatVec_otf proc~tmatvec_otf->proc~field_extd2d proc~elec_prop_d2m elec_prop_D2M proc~elec_prop_d2m->proc~prepare_fmm_ipd proc~prepare_polelec prepare_polelec proc~prepare_polelec->proc~elec_prop_m2d proc~prepare_polelec->proc~elec_prop_d2m proc~elec_prop_d2d elec_prop_D2D proc~prepare_polelec->proc~elec_prop_d2d proc~prepare_fixedelec prepare_fixedelec proc~prepare_fixedelec->proc~elec_prop_m2m proc~elec_prop_d2d->proc~prepare_fmm_ipd proc~polelec_geomgrad polelec_geomgrad proc~polelec_geomgrad->proc~prepare_polelec proc~ommp_get_polelec_energy ommp_get_polelec_energy proc~ommp_get_polelec_energy->proc~prepare_polelec proc~energy_mm_pol energy_MM_pol proc~ommp_get_polelec_energy->proc~energy_mm_pol proc~ommp_set_external_field ommp_set_external_field proc~ommp_set_external_field->proc~prepare_polelec proc~fixedelec_geomgrad fixedelec_geomgrad proc~fixedelec_geomgrad->proc~prepare_fixedelec proc~energy_mm_pol->proc~prepare_polelec proc~energy_mm_mm energy_MM_MM proc~energy_mm_mm->proc~prepare_fixedelec proc~ommp_polelec_geomgrad ommp_polelec_geomgrad proc~ommp_polelec_geomgrad->proc~polelec_geomgrad proc~c_ommp_get_polelec_energy C_ommp_get_polelec_energy proc~c_ommp_get_polelec_energy->proc~ommp_get_polelec_energy proc~ommp_get_full_ele_energy ommp_get_full_ele_energy proc~ommp_get_full_ele_energy->proc~ommp_get_polelec_energy proc~ommp_get_fixedelec_energy ommp_get_fixedelec_energy proc~ommp_get_full_ele_energy->proc~ommp_get_fixedelec_energy proc~ommp_set_external_field_nomm ommp_set_external_field_nomm proc~ommp_set_external_field_nomm->proc~ommp_set_external_field proc~ommp_full_geomgrad ommp_full_geomgrad proc~ommp_full_geomgrad->proc~polelec_geomgrad proc~ommp_full_geomgrad->proc~fixedelec_geomgrad proc~c_ommp_set_external_field_nomm C_ommp_set_external_field_nomm proc~c_ommp_set_external_field_nomm->proc~ommp_set_external_field proc~c_ommp_set_external_field C_ommp_set_external_field proc~c_ommp_set_external_field->proc~ommp_set_external_field proc~ommp_fixedelec_geomgrad ommp_fixedelec_geomgrad proc~ommp_fixedelec_geomgrad->proc~fixedelec_geomgrad proc~ommp_get_fixedelec_energy->proc~energy_mm_mm proc~c_ommp_polelec_geomgrad C_ommp_polelec_geomgrad proc~c_ommp_polelec_geomgrad->proc~ommp_polelec_geomgrad proc~c_ommp_full_geomgrad C_ommp_full_geomgrad proc~c_ommp_full_geomgrad->proc~ommp_full_geomgrad proc~c_ommp_get_full_ele_energy C_ommp_get_full_ele_energy proc~c_ommp_get_full_ele_energy->proc~ommp_get_full_ele_energy proc~c_ommp_fixedelec_geomgrad C_ommp_fixedelec_geomgrad proc~c_ommp_fixedelec_geomgrad->proc~ommp_fixedelec_geomgrad proc~ommp_get_full_energy ommp_get_full_energy proc~ommp_get_full_energy->proc~ommp_get_full_ele_energy proc~c_ommp_get_fixedelec_energy C_ommp_get_fixedelec_energy proc~c_ommp_get_fixedelec_energy->proc~ommp_get_fixedelec_energy proc~c_ommp_get_full_energy C_ommp_get_full_energy proc~c_ommp_get_full_energy->proc~ommp_get_full_energy

Contents

Source Code


Source Code

subroutine fmm_m2l_rotation_work(c, pm, pl, vscales, &
    & m2l_ztranslate_coef, alpha, src_m, beta, dst_l, work)
    ! Inputs
    integer, intent(in) :: pm, pl
    real(rp), intent(in) :: c(3), vscales((pm+pl+1)**2), &
        & m2l_ztranslate_coef(pm+1, pl+1, pl+1), alpha, src_m((pm+1)*(pm+1)), &
        & beta
    ! Output
    real(rp), intent(inout) :: dst_l((pl+1)*(pl+1))
    ! Temporary workspace
    real(rp), intent(out), target :: &
        & work(6*max(pm, pl)**2 + 19*max(pm, pl) + 8)
    ! Local variables
    real(rp) :: rho, ctheta, stheta, cphi, sphi
    integer :: m, n, p
    ! Pointers for temporary values of harmonics
    real(rp), pointer :: tmp_ml(:), tmp_ml2(:), vcos(:), vsin(:)
    ! Covert Cartesian coordinates into spherical
    call carttosph(c, rho, ctheta, stheta, cphi, sphi)
    ! If no need for rotations, just do translation along z
    if (abs(stheta) < eps_rp) then
        ! Workspace here is (pm+2)*(pm+1)
        call fmm_m2l_ztranslate_work(c(3), pm, pl, vscales, &
            & m2l_ztranslate_coef, alpha, src_m, beta, dst_l, work)
        return
    end if
    ! Prepare pointers
    p = max(pm, pl)
    m = (p+1)**2
    n = 4*m + 5*p ! 4*p*p + 13*p + 4
    tmp_ml(1:m) => work(n+1:n+m) ! 5*p*p + 15*p + 5
    n = n + m
    tmp_ml2(1:m) => work(n+1:n+m) ! 6*p*p + 17*p + 6
    n = n + m
    m = p + 1
    vcos => work(n+1:n+m) ! 6*p*p + 18*p + 7
    n = n + m
    vsin => work(n+1:n+m) ! 6*p*p + 19*p + 8
    ! Compute arrays of cos and sin that are needed for rotations of harmonics
    call trgev(cphi, sphi, p, vcos, vsin)
    ! Rotate around OZ axis (work array might appear in the future)
    call fmm_sph_rotate_oz_adj_work(pm, vcos, vsin, alpha, src_m, 0.0_rp, tmp_ml)
    ! Perform rotation in the OXZ plane, work size is 4*pm*pm+13*pm+4
    call fmm_sph_rotate_oxz_work(pm, ctheta, -stheta, 1.0_rp, tmp_ml, 0.0_rp, &
        & tmp_ml2, work)
    ! OZ translation, workspace here is (pm+2)*(pm+1)
    call fmm_m2l_ztranslate_work(rho, pm, pl, vscales, &
        & m2l_ztranslate_coef, 1.0_rp, tmp_ml2, 0.0_rp, tmp_ml, work)
    ! Backward rotation in the OXZ plane, work size is 4*pl*pl+13*pl+4
    call fmm_sph_rotate_oxz_work(pl, ctheta, stheta, 1.0_rp, tmp_ml, 0.0_rp, &
        & tmp_ml2, work)
    ! Backward rotation around OZ axis (work array might appear in the future)
    call fmm_sph_rotate_oz_work(pl, vcos, vsin, 1.0_rp, tmp_ml2, beta, dst_l)
end subroutine fmm_m2l_rotation_work