From 02608cf25a502f6777f33583bf00a4752d0319f4 Mon Sep 17 00:00:00 2001 From: Adam Parler Date: Wed, 29 Nov 2023 23:29:08 -0800 Subject: [PATCH] Added some scf files --- src/ddscf/comp4_ext.c | 32 +++++++++ src/ddscf/fock_2e_file.c | 142 +++++++++++++++++++++++++++++++++++++++ 2 files changed, 174 insertions(+) create mode 100644 src/ddscf/comp4_ext.c create mode 100644 src/ddscf/fock_2e_file.c diff --git a/src/ddscf/comp4_ext.c b/src/ddscf/comp4_ext.c new file mode 100644 index 0000000..b97284e --- /dev/null +++ b/src/ddscf/comp4_ext.c @@ -0,0 +1,32 @@ +#include +#include "bitops_decls.h" +#include "bitops_funcs.h" + +void comp4_extract(int* m, int i, double s, int nb_per_i) { + + int v; // Value after compression + +#if defined(CRAY) + int vv, vvv; +#endif + + int index, nbits; + double fast[] = {0.0, 1.0e-13, 1.0e-12, 1.0e-11, 1.0e-10, 1.0e-9, + 1.0e-8, 1.0e-7, 1.0e-6, 1.0e-5, 1.0e-4, 1.0e-3, 1.0e-2, + 1.0e-1, 1.0e0, 1.0e1}; + + v = 15; + index = (i - 1)/(2*nb_per_i) + 1; + nbits = 4*(i - (index-1)*(2*nb_per_i) - 1); +#if defined(CRAY) + vvv = shiftl(v, nbits); + vv = shiftr(iand(m(index), vvv), nbits); + v = iand(vv,15); +#else + v = iand(ishft(iand(m(index), ishft(v, nbits)), -nbits),15); +#endif + + s = fast(v); + // printf('%d -> %d %d %d %0.4f'); + +} diff --git a/src/ddscf/fock_2e_file.c b/src/ddscf/fock_2e_file.c new file mode 100644 index 0000000..fa0021c --- /dev/null +++ b/src/ddscf/fock_2e_file.c @@ -0,0 +1,142 @@ +#include "util.h" +#include "cscfps.h" +#include "cfock.h" +#include +#include + +void fock_2e_from_file(int geom, int basis, int nfock, int ablklen, + double jfac[nfock], double kfac[nfock], double tol2e, bool oskel, + double dij[nfock*ablklen], double dik[nfock*ablklen], double dli[nfock*ablklen], + double djk[nfock*ablklen], double dlj[nfock*ablklen], double dlk[nfock*ablklen], + double fij[nfock*ablklen], double fik[nfock*ablklen], double fli[nfock*ablklen], + double fjk[nfock*ablklen], double flj[nfock*ablklen], double flk[nfock*ablklen], + double tmp, int vg_dens[nfock], int vg_fock[nfock]) { + + //$Id$ + + /*Accumulate the contribution to the fock matrices from + integrals store in the integral file. Simply read thru + the file getting a range of indices, fetch the corresponding + density matrix blocks and then read the integrals in that + block. + */ + + double den_tol, denmax, dtol2e; + int ilo, jlo, klo, llo; + int ihi, jhi, khi, lhi; + int ijk_prev[3][2]; + int blklen; + + bool int2e_get_bf_range, int2e_file_read; + + if (oscfps) pstat_on(ps_fock_io); + + den_tol = fmax(tol2e*0.01, 1e-300); // To avoid a hard zero + + ijk_prev[0][0] = -1; + ijk_prev[1][0] = -1; + ijk_prev[2][0] = -1; + ijk_prev[0][1] = -1; + ijk_prev[1][1] = -1; + ijk_prev[2][1] = -1; + + blklen = nfock*ablklen; + dfill(blklen, 0.0e0, fij, 1); + dfill(blklen, 0.0e0, fik, 1); + dfill(blklen, 0.0e0, fli, 1); + dfill(blklen, 0.0e0, fjk, 1); + dfill(blklen, 0.0e0, flj, 1); + dfill(blklen, 0.0e0, flk, 1); + + // Loop over blocks of integral labels + + while (int2e_get_bf_range(ilo,ihi,jlo,jhi,klo,khi,llo,lhi)) { + // Get matrices for this block of labels + fock_init_cmul(ihi-ilo+1,jhi-jlo+1,lhi-llo+1); + fock_2e_cache_dens_fock( + ilo, jlo, klo, llo, + ihi, jhi, khi, lhi, + ijk_prev, + nfock, vg_dens, vg_fock, + jfac, kfac, + dij, dik, dli, djk, dlj, dlk, + fij, fik, fli, fjk, flj, flk, + tmp); + + fock_density_screen(nfock, + ilo, jlo, klo, llo, + ihi, jhi, khi, lhi, + ilo, jlo, klo, llo, + ihi, jhi, khi, lhi, + dij, dik, dli, djk, dlj, dlk, denmax) + + dtol2e = min(dentolmax, den_tol/max(1e-10,denmax), den_tol/max(1e-10,denmax**2)) + + call int2e_file_fock_block(nfock, dtol2e, + dij, dik, dli, djk, dlj, dlk, + fij, fik, fli, fjk, flj, flk) + + // Update F blocks + + call fock_upd_blk(nfock, vg_fock, + llo, lhi, ilo, ihi, kfac, fli, tmp) + call fock_upd_blk(nfock, vg_fock, + llo, lhi, jlo, jhi, kfac, flj, tmp) + call fock_upd_blk(nfock, vg_fock, + llo, lhi, klo, khi, jfac, flk, tmp) + } + + if (ijk_prev[0][0]) != -1) { + fock_upd_blk(nfock, vg_fock, + ijk_prev[0][0]), ijk_prev[0][1]), + ijk_prev(2,1), ijk_prev(2,2), + jfac, fij, tmp) + fock_upd_blk(nfock, vg_fock, + ijk_prev(2,1), ijk_prev(2,2), + ijk_prev(3,1), ijk_prev(3,2), + kfac, fjk, tmp ) + fock_upd_blk( nfock, vg_fock, + ijk_prev(1,1), ijk_prev(1,2), + ijk_prev(3,1), ijk_prev(3,2), + kfac, fik, tmp ) + } + + if (oscfps) pstat_off(ps_fock_io); + +} + +void fock_2e_rep_from_file(int geom, int basis, int nfock, int nbf, + double jfac[nfock], double kfac[nfock], double tol2e, bool oskel, + double dens[nfock][nbf*nbf], fock[nfock][nbf*nbf]) { + + double den_tol, denmax; + int ilo, jlo, klo, llo; + int ihi, jhi, khi, lhi, i, j; + + bool int2e_get_bf_range, int2e_file_read; + int idamax; + + if (oscfps) pstat_on(ps_fock_io); + + denmax = 0.0; + for (i = 0; i < nfock; i++) { + j = idamax(nbf*nbf, dens[i][0], nfock); + denmax = max(denmax, abs(dens[i][j]); + } + // return if DM is null (e.g imaginary part of RTTDFT DM at t=0) + if (denmax < 1e-12) return; + den_tol = min(dentolmax,tol2e/denmax,tol2e/denmax**2) // Threshold to screen integs only + + if (ga_nodeid() == 0 && util_print('fockfile',print_debug)) { + printf("fockfile: tols %d %d %d %d", tol2e, dentolmax, denmax, den_tol); + } + + fock_init_cmul(nbf,nbf,nbf) // lookup table for f build + + // Loop over blocks of integral labels +} + + + + +}