################################### lib_bound: bound state calculations ################################### ******************************** Potential energy matrix elements ******************************** Potential matrix elements: atom + open-shell diatom, body fixed =============================================================== .. function:: qmsV = qms_V_bf_DLMK2(sjmkf, sjmki, sLMK, qnames) Compute matrix elements of the complex conjugate of a Wigner rotation matrix in a symmetric top basis. .. math:: \begin{align} V_{fi}(L_n,M_n,K_n) & = \langle j_f m_f k_f | D^{(L)\ast}_{MK}|j_i m_i k_i\rangle\\ &= \frac{[j_f,j_i]^{1/2}}{8\pi^2} \int_0^{2\pi} d\alpha \int_0^\pi \sin\beta d\beta \int_0^{2\pi} d\gamma\\ &~~~~\mbox{}\times D^{(j_f)}_{m_fk_f}(\alpha,\beta,\gamma) D^{(L_n)\ast}_{M_nK_n}(\alpha,\beta,\gamma) D^{(j_i)\ast}_{m_ik_i}(\alpha,\beta,\gamma)\,\\ &= [j_f,j_i]^{1/2} (-1)^{m_f-k_f} \begin{pmatrix} j_f & L_n & j_i\\ -m_f & M_n & m_i \end{pmatrix} \begin{pmatrix} j_f & L_n & j_i\\ -k_f & K_n & k_i \end{pmatrix}, \end{align} .. where :math:`[j_f,f_i] = (2j_f+1)(2j_i+1)` and :math:`D^{(L)}_{MK}` are Wigner D-matrices. This function implements Eqs. (11) and (19) in a study of He+HF :sup:`+` :cite:`dhont:04`. It can also be used for Eq. (20) in a study of He+CO :sub:`2` :cite:`selim:23`. Note that a factor :math:`\sqrt{(2j'+1)(2j+1)}` is missing in that equation. **Input arguments** :sjmkf: A ``q_set`` data structure (see ``lib_qmat``) containing the set of bra quantum numbers: :math:`\{j_f m_f k_f; f=1,2,\ldots\}` Or: A Scilab ``matrix`` with two columns, with ``sjmkf(f,:) = [j_f m_f k_f]`` :sjmki: A ``q_set`` data structure (see ``lib_qmat``) containing the set of ket quantum numbers :math:`\{j_i m_i k_i; i=1,2,\ldots\}` Or: A Scilab ``matrix`` with two columns, with ``sjmki(i,:) = [j_i m_i k_i]`` :sLMK: A ``q_set`` data structure (see ``lib_qmat``) containing :math:`LMK` of the Wigner *D*-matrix element :math:`D^{(L)\ast}_{MK}`. Or: A Scilab :math:`n \times 3` ``matrix`` ``Q`` with elements :math:`Q(i,:)=[L_i M_i K_i]` :qnames: A ``vector`` of ``strings`` with the names of the quantum numbers. By default ``qnames = ['j' 'm' 'k' 'L' 'M' 'K']`` **Output argument** :qmsV: A ``q_mats`` data structure (see ``lib_qmat``) containing a list of matrices with Gaunt coefficients for all values of :math:`L` .. code-block:: Scilab qmsV.s1 = sjmkf qmsV.s2 = sjmki qmsV.sa = sLMK qmsV.As(i) = V_n // list of matrices V_n(f, i) = **Example:** Here is an example of how to use this function: .. code-block:: Scilab // Example usage of the function qms_V_bf_PL2 lib_load(['lib_bound' 'lib_special' 'lib_tenop']) jmk = [0 0 0 1 0 0 1 0 1 2 0 0 2 0 1 2 0 2] LMKs = [0 0 0 1 0 0 1 0 1 2 0 0 2 0 1 2 0 2] qmsV1 = qms_V_bf_DLMK2(jmk, jmk, LMKs, ['j' 'm' 'k' 'L' 'M' 'K']) // |J(3), K(3)>A(6;6)<(same)| [L(5); 5] // Or, create q_sets first s1 = q_set(jmk, ['j' 'm' 'k']) sL = q_set(LMKs, ['L' 'M' 'K']) // LMKs has 3 columns qmsV2 = qms_V_bf_DLMK2(s1, s1, sL, [s1.qnames sL.qnames]) if qmsV1==qmsV2 mprintf('The same!\n') else mprintf('qmsV1 <> qmsV2\n') end Gaunt coefficients ================== .. function:: qmsV = qms_V_bf_PL2(sjkf, sjki, sL, qnames) This function computes matrices with Gaunt coefficients: .. math:: \begin{align} V_{fi}(L) & = \langle j_f k_f | P_L|j_i k_i\rangle\\ &=\int_0^\pi \! \int_0^{2\pi} Y^\ast_{j_fk_f}(\theta,\phi) P_L(\cos\theta) Y_{j_ik_i}(\theta,\phi)\, d\phi \, \sin\theta d\theta\\ &= [(2j_f+1)(2j_i+1)]^{1/2} (-1)^k \begin{pmatrix} j_f & L & j_i\\ 0 & 0 & 0 \end{pmatrix} \begin{pmatrix} j_f & L & j_i\\ -k_f & 0 & k_i \end{pmatrix}, \end{align} .. where :math:`Y_{jk}` are spherical harmonics with the Condon and Shortley phase convention and :math:`P_L` are Legendre polynomials. **Input arguments** :sjkf: A ``q_set`` data structure (see ``lib_qmat``) containing the set of bra quantum numbers: :math:`\{j_f k_f; f=1,2,\ldots\}` Or: A Scilab ``matrix`` with two columns, with ``sjkf(f,:) = [j_f k_f]`` :sjki: A ``q_set`` data structure (see ``lib_qmat``) containing the set of ket quantum numbers :math:`\{j_i k_i; i=1,2,\ldots\}` Or: A Scilab ``matrix`` with two columns, with ``sjki(i,:) = [j_i k_i]`` :sL: A ``q_set`` data structure (see ``lib_qmat``) containing the orders :math:`L` or the Legendre polynomials :math:`P_L` Or: A Scilab ``vector`` with the values of :math:`L` :qnames: A ``vector`` of ``strings`` with the names of the quantum numbers. By default ``qnames = ['j' 'k' 'L']`` **Output argument** :qmsV: A ``q_mats`` data structure (see ``lib_qmat``) containing a list of matrices with Gaunt coefficients for all values of :math:`L` .. code-block:: Scilab qmsV.s1 = sjkf qmsV.s2 = sjki qmsV.sa = sL qmsV.As(i) = VL // list of Gaunt coefficient matrices for: L = sL(i, 'L') // VL(f, i) = **Example:** Here is an example of how to use this function: .. code-block:: Scilab // Example usage of the function qms_V_bf_PL2 lib_load(['lib_bound' 'lib_special' 'lib_tenop']) jk = [0 0 1 0 1 1 2 0 2 1 2 2] Ls = 0:4 qmsV1 = qms_V_bf_PL2(jk, jk, Ls, ['j' 'k' 'L']) // |J(3), K(3)>A(6;6)<(same)| [L(5); 5] // Or, create q_sets first s1 = q_set(jk, ['j' 'k']) sL = q_set(Ls(:), 'L') // q_set requires a column vector here qmsV2 = qms_V_bf_PL2(s1, s1, sL, [s1.qnames 'L']) if qmsV1==qmsV2 mprintf('The same!\n') else mprintf('qmsV1 <> qmsV2\n') end .. bibliography::