1
0
mirror of https://gitlab.com/scemama/qp_plugins_scemama.git synced 2024-11-07 22:53:42 +01:00
qp_plugins_scemama/devel/trexio/export_trexio.irp.f
2021-10-12 16:00:27 +02:00

276 lines
6.6 KiB
Fortran

program export_trexio
use trexio
implicit none
BEGIN_DOC
! Exports the wave function in TREXIO format
END_DOC
integer(8) :: f ! TREXIO file handle
integer :: rc
print *, 'TREXIO file : '//trim(trexio_filename)
print *, ''
call system('rm -rf '//trim(trexio_filename)//'.bak')
call system('mv '//trim(trexio_filename)//' '//trim(trexio_filename)//'.bak')
! call system('mv '//trim(trexio_filename)//' '//trim(trexio_filename)''//.bak)
if (backend == 0) then
f = trexio_open(trexio_filename, 'w', TREXIO_HDF5, rc)
else if (backend == 1) then
f = trexio_open(trexio_filename, 'w', TREXIO_TEXT, rc)
endif
if (f == 0) then
print *, 'Unable to open TREXIO file for writing'
print *, 'rc = ', rc
stop -1
endif
! ------------------------------------------------------------------------------
! Electrons
! ---------
rc = trexio_write_electron_up_num(f, elec_alpha_num)
call check_success(rc)
rc = trexio_write_electron_dn_num(f, elec_beta_num)
call check_success(rc)
! Nuclei
! ------
rc = trexio_write_nucleus_num(f, nucl_num)
call check_success(rc)
rc = trexio_write_nucleus_charge(f, nucl_charge)
call check_success(rc)
rc = trexio_write_nucleus_coord(f, nucl_coord_transp)
call check_success(rc)
rc = trexio_write_nucleus_label(f, nucl_label, 32)
call check_success(rc)
! Pseudo-potentials
! -----------------
double precision, allocatable :: tmp_double(:,:)
integer, allocatable :: tmp_int(:,:)
! rc = trexio_write_ecp_lmax_plus_1(f, pseudo_lmax+1)
! call check_success(rc)
!
! rc = trexio_write_ecp_z_core(f, nucl_charge_remove)
! call check_success(rc)
!
! rc = trexio_write_ecp_local_num_n_max(f, pseudo_klocmax)
! call check_success(rc)
!
! rc = trexio_write_ecp_local_power(f, pseudo_n_k_transp)
! call check_success(rc)
!
! rc = trexio_write_ecp_local_exponent(f, pseudo_dz_k_transp)
! call check_success(rc)
!
! rc = trexio_write_ecp_local_coef(f, pseudo_v_k_transp)
! call check_success(rc)
!
! rc = trexio_write_ecp_non_local_num_n_max(f, pseudo_kmax)
! call check_success(rc)
!
! rc = trexio_write_ecp_non_local_power(f, pseudo_n_kl_transp)
! call check_success(rc)
!
! rc = trexio_write_ecp_non_local_exponent(f, pseudo_dz_kl_transp)
! call check_success(rc)
!
! rc = trexio_write_ecp_non_local_coef(f, pseudo_v_kl_transp)
! call check_success(rc)
! Basis
! -----
rc = trexio_write_basis_type(f, 'Gaussian', len('Gaussian'))
call check_success(rc)
rc = trexio_write_basis_num(f, shell_num)
call check_success(rc)
rc = trexio_write_basis_nucleus_shell_num(f, nucleus_shell_num)
call check_success(rc)
rc = trexio_write_basis_nucleus_index(f, basis_nucleus_index)
call check_success(rc)
rc = trexio_write_basis_shell_ang_mom(f, shell_ang_mom)
call check_success(rc)
rc = trexio_write_basis_prim_num(f, prim_num)
call check_success(rc)
rc = trexio_write_basis_shell_prim_num(f, shell_prim_num)
call check_success(rc)
double precision, allocatable :: factor(:)
allocate(factor(shell_num))
if (ao_normalized) then
factor(1:shell_num) = shell_normalization_factor(1:shell_num)
else
factor(1:shell_num) = 1.d0
endif
rc = trexio_write_basis_shell_factor(f, factor)
call check_success(rc)
deallocate(factor)
rc = trexio_write_basis_shell_prim_index(f, shell_prim_index)
call check_success(rc)
rc = trexio_write_basis_exponent(f, prim_expo)
call check_success(rc)
rc = trexio_write_basis_coefficient(f, prim_coef)
call check_success(rc)
allocate(factor(prim_num))
if (primitives_normalized) then
factor(1:prim_num) = prim_normalization_factor(1:prim_num)
else
factor(1:prim_num) = 1.d0
endif
rc = trexio_write_basis_prim_factor(f, factor)
call check_success(rc)
deallocate(factor)
! Atomic orbitals
! ---------------
rc = trexio_write_ao_num(f, ao_num)
call check_success(rc)
rc = trexio_write_ao_cartesian(f, 1)
call check_success(rc)
rc = trexio_write_ao_shell(f, ao_shell)
call check_success(rc)
integer :: i, pow0(3), powA(3), j, k, l, nz
double precision :: normA, norm0, C_A(3), overlap_x, overlap_z, overlap_y, c
nz=100
C_A(1) = 0.d0
C_A(2) = 0.d0
C_A(3) = 0.d0
allocate(factor(ao_num))
print *, ao_first_of_shell
if (ao_normalized) then
do i=1,ao_num
l = ao_first_of_shell(ao_shell(i))
factor(i) = (ao_coef_normalized(i,1)+tiny(1.d0))/(ao_coef_normalized(l,1)+tiny(1.d0))
enddo
else
factor(:) = 1.d0
endif
rc = trexio_write_ao_normalization(f, factor)
call check_success(rc)
deallocate(factor)
! One-e AO integrals
! ------------------
rc = trexio_write_ao_1e_int_overlap(f,ao_overlap)
call check_success(rc)
rc = trexio_write_ao_1e_int_kinetic(f,ao_kinetic_integrals)
call check_success(rc)
rc = trexio_write_ao_1e_int_potential_n_e(f,ao_integrals_n_e)
call check_success(rc)
if (do_pseudo) then
rc = trexio_write_ao_1e_int_ecp_local(f,ao_pseudo_integrals_local)
call check_success(rc)
rc = trexio_write_ao_1e_int_ecp_non_local(f,ao_pseudo_integrals_non_local)
call check_success(rc)
endif
rc = trexio_write_ao_1e_int_core_hamiltonian(f,ao_one_e_integrals)
call check_success(rc)
! Molecular orbitals
! ------------------
! rc = trexio_write_mo_type(f, mo_label)
! call check_success(rc)
rc = trexio_write_mo_num(f, mo_num)
call check_success(rc)
rc = trexio_write_mo_coefficient(f, mo_coef)
call check_success(rc)
! One-e MO integrals
! ------------------
rc = trexio_write_mo_1e_int_kinetic(f,mo_kinetic_integrals)
call check_success(rc)
rc = trexio_write_mo_1e_int_potential_n_e(f,mo_integrals_n_e)
call check_success(rc)
if (do_pseudo) then
rc = trexio_write_mo_1e_int_ecp_local(f,mo_pseudo_integrals_local)
call check_success(rc)
rc = trexio_write_mo_1e_int_ecp_non_local(f,mo_pseudo_integrals_non_local)
call check_success(rc)
endif
!
rc = trexio_write_mo_1e_int_core_hamiltonian(f,one_e_dm_mo)
call check_success(rc)
! RDM
! ----
! rc = trexio_write_rdm_one_e(f,one_e_dm_mo)
! call check_success(rc)
!
! rc = trexio_write_rdm_one_e_up(f,one_e_dm_mo_alpha_average)
! call check_success(rc)
!
! rc = trexio_write_rdm_one_e_dn(f,one_e_dm_mo_beta_average)
! call check_success(rc)
! ------------------------------------------------------------------------------
rc = trexio_close(f)
call check_success(rc)
end
subroutine check_success(rc)
use trexio
implicit none
integer, intent(in) :: rc
character*(128) :: str
if (rc /= TREXIO_SUCCESS) then
call trexio_string_of_error(rc,str)
print *, 'TREXIO Error: ' //trim(str)
stop -1
endif
end
! -*- mode: f90 -*-