# Path: orbkit/read/high_level.py
# def main_read(fname, all_mo=False, spin=None, itype='auto', check_norm=False, **kwargs):
#   '''
#   This is the high-lever interface for the
#   orbkit reading routines.
#   
#   **Parameters:**
#   
#   fname: str, file descriptor
#     Specifies the filename for the input file.
#     fname can also be used with a file descriptor instad of a filename.
#   all_mo : bool, optional
#     If True, all molecular orbitals are returned.
#   spin : {None, 'alpha', or 'beta'}, optional
#     If not None, returns exclusively 'alpha' or 'beta' molecular orbitals.
#   itype : str, optional
#     Can be used to manually specify the input filetype.
#   check_norm : bool, optional
#     If True, ORBKIT verifies that molecular orbitals are orthonormal.
#     
#   **Note:**
#   
#     All additional keyword arguments are forwarded to the reading functions.
#     
#   **Returns:**
#   
#     qc (class QCinfo) with attributes geo_spec, geo_info, ao_spec, mo_spec, etot :
#     See :ref:`Central Variables` for details.
#   '''
# 
#   if isinstance(fname, str):
#     filename = fname
#   else:
#     filename = fname.name
# 
#   if itype == 'auto':
#     itype = find_itype(fname)
#  
#   display('Loading data from {0} type file {1}\n'.format(itype, filename))
# 
#   qc = readers[itype](fname, all_mo=all_mo, spin=spin, **kwargs)
# 
#   if check_norm:
#     deviation = check_mo_norm(qc)
#     if deviation >= 1e-5:
#       raise ValueError('Bad molecular orbital norm: {0:%4e}'.format(deviation))
# 
#   return qc
# 
# Path: orbkit/analytical_integrals.py
# def get_ao_overlap(coord_a, coord_b, ao_spec, lxlylz_b=None,
#                    drv=None):
# def cartesian2spherical_aoom(ao_overlap_matrix,ao_spec):
# def get_mo_overlap(mo_a,mo_b,ao_overlap_matrix):
# def initializer(gargs):
# def get_slice(x):  
# def get_mo_overlap_matrix(mo_a,mo_b,ao_overlap_matrix,numproc=1):
# def get_moom_atoms(atoms,qc,mo_a,mo_b,ao_overlap_matrix,numproc=1):
# def check_mo_norm(qc):
# def get_dipole_moment(qc,component=['x','y','z']):
# def get_ao_dipole_matrix(qc,component='x'):
# def get_nuclear_dipole_moment(qc,component='x'):
# def get_atom2mo(qc):
# def get_lc(atoms,atom2mo,strict=False):


from orbkit.read import main_read
from orbkit.tools import *
from orbkit.analytical_integrals import get_ao_overlap,get_mo_overlap,print2D
This file is part of orbkit. See the main program or documentation 
for information on the license.

It shows two examples of how to use orbkit for evaluating 
analytical integrals.

First Part:

  Compute the electron number analytically using the standard density formula.

Second Part:

  Compute the dipole moments analytically. This expectation value is more complex.
  As an explanation, consider the following expectation value between two primitive 
  $p_z$ orbitals:

  .. math::

    <\phi_k|z|\phi_l> = \int dxdydz (z - Z_k)^{k_z} e^{\alpha_k r_k^2} (z - Z_l)^{l_z+1} e^{\alpha_l r_l^2} 
                        + Z_l \int dxdydz (z - Z_k)^{k_z} e^{\alpha_k r_k^2} (z - Z_l)^{l_z} e^{\alpha_l r_l^2} 
                      = ao_part_1 + ao_part_2

with :math:`r_l = \sqrt{(x-X_l)^2 + (y-Y_l)^2 + (z - Z_l)^2)}`
'''

in_fid = 'h2o.molden'
# Read the input file
qc = main_read(in_fid,itype='molden',all_mo=False)

# Compute atomic orbital overlap matrix 

# Path: orbkit/read/high_level.py
# def main_read(fname, all_mo=False, spin=None, itype='auto', check_norm=False, **kwargs):
#   '''
#   This is the high-lever interface for the
#   orbkit reading routines.
#   
#   **Parameters:**
#   
#   fname: str, file descriptor
#     Specifies the filename for the input file.
#     fname can also be used with a file descriptor instad of a filename.
#   all_mo : bool, optional
#     If True, all molecular orbitals are returned.
#   spin : {None, 'alpha', or 'beta'}, optional
#     If not None, returns exclusively 'alpha' or 'beta' molecular orbitals.
#   itype : str, optional
#     Can be used to manually specify the input filetype.
#   check_norm : bool, optional
#     If True, ORBKIT verifies that molecular orbitals are orthonormal.
#     
#   **Note:**
#   
#     All additional keyword arguments are forwarded to the reading functions.
#     
#   **Returns:**
#   
#     qc (class QCinfo) with attributes geo_spec, geo_info, ao_spec, mo_spec, etot :
#     See :ref:`Central Variables` for details.
#   '''
# 
#   if isinstance(fname, str):
#     filename = fname
#   else:
#     filename = fname.name
# 
#   if itype == 'auto':
#     itype = find_itype(fname)
#  
#   display('Loading data from {0} type file {1}\n'.format(itype, filename))
# 
#   qc = readers[itype](fname, all_mo=all_mo, spin=spin, **kwargs)
# 
#   if check_norm:
#     deviation = check_mo_norm(qc)
#     if deviation >= 1e-5:
#       raise ValueError('Bad molecular orbital norm: {0:%4e}'.format(deviation))
# 
#   return qc
# 
# Path: orbkit/analytical_integrals.py
# def get_ao_overlap(coord_a, coord_b, ao_spec, lxlylz_b=None,
#                    drv=None):
# def cartesian2spherical_aoom(ao_overlap_matrix,ao_spec):
# def get_mo_overlap(mo_a,mo_b,ao_overlap_matrix):
# def initializer(gargs):
# def get_slice(x):  
# def get_mo_overlap_matrix(mo_a,mo_b,ao_overlap_matrix,numproc=1):
# def get_moom_atoms(atoms,qc,mo_a,mo_b,ao_overlap_matrix,numproc=1):
# def check_mo_norm(qc):
# def get_dipole_moment(qc,component=['x','y','z']):
# def get_ao_dipole_matrix(qc,component='x'):
# def get_nuclear_dipole_moment(qc,component='x'):
# def get_atom2mo(qc):
# def get_lc(atoms,atom2mo,strict=False):


# Path: examples/basic_examples/orbkit_non-grid-based.py
from orbkit.read import main_read
from orbkit.tools import *
from orbkit.analytical_integrals import get_ao_overlap,get_mo_overlap,print2D

This file is part of orbkit. See the main program or documentation 
for information on the license.

It shows two examples of how to use orbkit for evaluating 
analytical integrals.

First Part:

  Compute the electron number analytically using the standard density formula.

Second Part:

  Compute the dipole moments analytically. This expectation value is more complex.
  As an explanation, consider the following expectation value between two primitive 
  $p_z$ orbitals:

  .. math::

    <\phi_k|z|\phi_l> = \int dxdydz (z - Z_k)^{k_z} e^{\alpha_k r_k^2} (z - Z_l)^{l_z+1} e^{\alpha_l r_l^2} 
                        + Z_l \int dxdydz (z - Z_k)^{k_z} e^{\alpha_k r_k^2} (z - Z_l)^{l_z} e^{\alpha_l r_l^2} 
                      = ao_part_1 + ao_part_2

with :math:`r_l = \sqrt{(x-X_l)^2 + (y-Y_l)^2 + (z - Z_l)^2)}`
'''

in_fid = 'h2o.molden'
# Read the input file
qc = main_read(in_fid,itype='molden',all_mo=False)

# Compute atomic orbital overlap matrix 

ao_overlap_matrix = get_ao_overlap(qc.geo_spec,qc.geo_spec,qc.ao_spec)