diff --git a/Makefile.am b/Makefile.am index 8fae847..1137766 100644 --- a/Makefile.am +++ b/Makefile.am @@ -47,9 +47,9 @@ pkgconfig_DATA = pkgconfig/qmckl.pc qmckl_h = include/qmckl.h include_HEADERS = $(qmckl_h) -test_qmckl_f = tests/qmckl_f.f90 +test_qmckl_f = tests/qmckl_f.F90 test_qmckl_fo = tests/qmckl_f.o -src_qmckl_f = src/qmckl_f.f90 +src_qmckl_f = src/qmckl_f.F90 src_qmckl_fo = src/qmckl_f.o header_tests = tests/chbrclf.h tests/n2.h @@ -139,7 +139,7 @@ cat_h_verbose_0 = @echo " HEADER $@"; ## Rules ## ===== -SUFFIXES = .f90 .h .org .c _f.f90 _func.h _type.h _private_func.h _private_type.h +SUFFIXES = .F90 .h .org .c _f.F90 _func.h _type.h _private_func.h _private_type.h $(test_qmckl_f): $(src_qmckl_f) cp $(src_qmckl_f) $(test_qmckl_f) diff --git a/configure.ac b/configure.ac index 1f58783..6654fd7 100644 --- a/configure.ac +++ b/configure.ac @@ -71,13 +71,16 @@ AC_LANG(C) # Checks for programs. AC_PROG_CC +AC_PROG_F77 + # Make sure the c compiler supports C99 m4_version_prereq([2.70],[], [AC_PROG_CC_C99]) AS_IF([test "$ac_cv_prog_cc_c99" = "no"], [AC_MSG_ERROR([The compiler does not support C99])]) AC_PROG_CC_C_O AC_PROG_FC AC_PROG_FC_C_O -AC_FC_SRCEXT([f90]) +AC_FC_PP_DEFINE +AC_FC_SRCEXT([F90]) AC_FC_FREEFORM LT_INIT AC_PROG_INSTALL @@ -192,6 +195,10 @@ esac # Options. +AC_ARG_ENABLE(hpc, [AS_HELP_STRING([--enable-hpc],[Use HPC-optimized functions])], HAVE_HPC=$enableval, HAVE_HPC=no) +AS_IF([test "$HAVE_HPC" = "yes"], [ + AC_DEFINE([HAVE_HPC], [1], [If defined, activate HPC routines]) +]) AC_ARG_ENABLE(debug, [AS_HELP_STRING([--enable-debug],[compile for debugging])], ok=$enableval, ok=no) if test "$ok" = "yes"; then @@ -336,6 +343,7 @@ FCLAGS..........: ${FCFLAGS} LDFLAGS:........: ${LDFLAGS} LIBS............: ${LIBS} USE CHAMELEON...: ${with_chameleon} +HPC version.....: ${HAVE_HPC} Package features: ${ARGS} diff --git a/org/qmckl_ao.org b/org/qmckl_ao.org index 28a831b..d40e478 100644 --- a/org/qmckl_ao.org +++ b/org/qmckl_ao.org @@ -58,6 +58,12 @@ gradients and Laplacian of the atomic basis functions. #include #+end_src + #+begin_src f90 :tangle (eval f) :noweb yes +#ifdef HAVE_CONFIG_H +#include "config.h" +#endif + #+end_src + #+begin_src c :tangle (eval c_test) :noweb yes #include "qmckl.h" #include "assert.h" @@ -2460,11 +2466,11 @@ for (int64_t i=0 ; i < ao_num ; ++i) { |----------------------+-----------------------------------+----------------------------------------------------------------------------------------------| | Variable | Type | Description | |----------------------+-----------------------------------+----------------------------------------------------------------------------------------------| - | ~primitive_vgl~ | ~double[5][point_num][prim_num]~ | Value, gradients, Laplacian of the primitives at current positions | + | ~primitive_vgl~ | ~double[point_num][5][prim_num]~ | Value, gradients, Laplacian of the primitives at current positions | | ~primitive_vgl_date~ | ~uint64_t~ | Last modification date of Value, gradients, Laplacian of the primitives at current positions | - | ~shell_vgl~ | ~double[5][point_num][shell_num]~ | Value, gradients, Laplacian of the primitives at current positions | + | ~shell_vgl~ | ~double[point_num][5][shell_num]~ | Value, gradients, Laplacian of the primitives at current positions | | ~shell_vgl_date~ | ~uint64_t~ | Last modification date of Value, gradients, Laplacian of the AOs at current positions | - | ~ao_vgl~ | ~double[5][point_num][ao_num]~ | Value, gradients, Laplacian of the primitives at current positions | + | ~ao_vgl~ | ~double[point_num][5][ao_num]~ | Value, gradients, Laplacian of the primitives at current positions | | ~ao_vgl_date~ | ~uint64_t~ | Last modification date of Value, gradients, Laplacian of the AOs at current positions | @@ -3055,7 +3061,7 @@ assert(0 == test_qmckl_ao_gaussian_vgl(context)); | ~coord~ | ~double[3][point_num]~ | in | Coordinates | | ~nucl_coord~ | ~double[3][nucl_num]~ | in | Nuclear coordinates | | ~expo~ | ~double[prim_num]~ | in | Exponents of the primitives | - | ~primitive_vgl~ | ~double[5][point_num][prim_num]~ | out | Value, gradients and Laplacian of the primitives | + | ~primitive_vgl~ | ~double[point_num][5][prim_num]~ | out | Value, gradients and Laplacian of the primitives | #+CALL: generate_c_header(table=qmckl_ao_basis_primitive_gaussian_vgl_args,rettyp=get_value("CRetType"),fname="qmckl_compute_ao_basis_primitive_gaussian_vgl")) @@ -3091,7 +3097,7 @@ integer function qmckl_compute_ao_basis_primitive_gaussian_vgl_f( & double precision , intent(in) :: coord(point_num,3) double precision , intent(in) :: nucl_coord(nucl_num,3) double precision , intent(in) :: expo(prim_num) - double precision , intent(out) :: primitive_vgl(prim_num,point_num,5) + double precision , intent(out) :: primitive_vgl(prim_num,5,point_num) integer*8 :: inucl, iprim, ipoint double precision :: x, y, z, two_a, ar2, r2, v, cutoff @@ -3116,11 +3122,11 @@ integer function qmckl_compute_ao_basis_primitive_gaussian_vgl_f( & v = dexp(-ar2) two_a = -2.d0 * expo(iprim) * v - primitive_vgl(iprim, ipoint, 1) = v - primitive_vgl(iprim, ipoint, 2) = two_a * x - primitive_vgl(iprim, ipoint, 3) = two_a * y - primitive_vgl(iprim, ipoint, 4) = two_a * z - primitive_vgl(iprim, ipoint, 5) = two_a * (3.d0 - 2.d0*ar2) + primitive_vgl(iprim, 1, ipoint) = v + primitive_vgl(iprim, 2, ipoint) = two_a * x + primitive_vgl(iprim, 3, ipoint) = two_a * y + primitive_vgl(iprim, 4, ipoint) = two_a * z + primitive_vgl(iprim, 5, ipoint) = two_a * (3.d0 - 2.d0*ar2) end do end do @@ -3156,7 +3162,7 @@ end function qmckl_compute_ao_basis_primitive_gaussian_vgl_f real (c_double ) , intent(in) :: coord(point_num,3) real (c_double ) , intent(in) :: nucl_coord(nucl_num,3) real (c_double ) , intent(in) :: expo(prim_num) - real (c_double ) , intent(out) :: primitive_vgl(prim_num,point_num,5) + real (c_double ) , intent(out) :: primitive_vgl(prim_num,5,point_num) integer(c_int32_t), external :: qmckl_compute_ao_basis_primitive_gaussian_vgl_f info = qmckl_compute_ao_basis_primitive_gaussian_vgl_f & @@ -3317,17 +3323,17 @@ print ( "[7][4][26] : %e"% lf(a,x,y)) - double prim_vgl[5][elec_num*walk_num][prim_num]; + double prim_vgl[elec_num*walk_num][5][prim_num]; rc = qmckl_get_ao_basis_primitive_vgl(context, &(prim_vgl[0][0][0]), (int64_t) 5*elec_num*walk_num*prim_num ); assert (rc == QMCKL_SUCCESS); - assert( fabs(prim_vgl[0][26][7] - ( 1.0501570432064878E-003)) < 1.e-14 ); - assert( fabs(prim_vgl[1][26][7] - (-7.5014974095310560E-004)) < 1.e-14 ); - assert( fabs(prim_vgl[2][26][7] - (-3.8250692897610380E-003)) < 1.e-14 ); - assert( fabs(prim_vgl[3][26][7] - ( 3.4950559194080275E-003)) < 1.e-14 ); - assert( fabs(prim_vgl[4][26][7] - ( 2.0392163767356572E-002)) < 1.e-14 ); + assert( fabs(prim_vgl[26][0][7] - ( 1.0501570432064878E-003)) < 1.e-14 ); + assert( fabs(prim_vgl[26][1][7] - (-7.5014974095310560E-004)) < 1.e-14 ); + assert( fabs(prim_vgl[26][2][7] - (-3.8250692897610380E-003)) < 1.e-14 ); + assert( fabs(prim_vgl[26][3][7] - ( 3.4950559194080275E-003)) < 1.e-14 ); + assert( fabs(prim_vgl[26][4][7] - ( 2.0392163767356572E-002)) < 1.e-14 ); } @@ -3387,7 +3393,7 @@ for (j=0 ; j 0) then + pows(-2:0,1:3) = 1.d0 + do i=1,lmax + pows(i,1) = pows(i-1,1) * Y(1) + pows(i,2) = pows(i-1,2) * Y(2) + pows(i,3) = pows(i-1,3) * Y(3) + end do + + l (1:3,1:4) = 0 + VGL(1:4,1:5) = 0.d0 + + VGL(1 ,1 ) = 1.d0 + VGL(2:4,1:5) = 0.d0 + + l (1,2) = 1 + VGL(2,1) = pows(1,1) + VGL(2,2) = 1.d0 + + l (2,3) = 1 + VGL(3,1) = pows(1,2) + VGL(3,3) = 1.d0 + + l (3,4) = 1 + VGL(4,1) = pows(1,3) + VGL(4,4) = 1.d0 + + n=4 + endif + + ! l>=2 + dd = 2.d0 + do d=2,lmax + da = dd + do a=d,0,-1 + db = dd-da + do b=d-a,0,-1 + c = d - a - b + dc = dd - da - db + n = n+1 + + l(1,n) = a + l(2,n) = b + l(3,n) = c + + xy = pows(a,1) * pows(b,2) + yz = pows(b,2) * pows(c,3) + xz = pows(a,1) * pows(c,3) + + VGL(n,1) = xy * pows(c,3) + + xy = dc * xy + xz = db * xz + yz = da * yz + + VGL(n,2) = pows(a-1,1) * yz + VGL(n,3) = pows(b-1,2) * xz + VGL(n,4) = pows(c-1,3) * xy + + VGL(n,5) = & + (da-1.d0) * pows(a-2,1) * yz + & + (db-1.d0) * pows(b-2,2) * xz + & + (dc-1.d0) * pows(c-2,3) * xy + + db = db - 1.d0 + end do + da = da - 1.d0 + end do + dd = dd + 1.d0 + end do + + info = QMCKL_SUCCESS + +end function qmckl_ao_polynomial_transp_vgl_f + #+end_src + + #+CALL: generate_c_interface(table=qmckl_ao_polynomial_vgl_args,rettyp=get_value("CRetType"),fname=get_value("Name")) + + #+RESULTS: + #+begin_src f90 :tangle (eval f) :comments org :exports none + integer(c_int32_t) function qmckl_ao_polynomial_transp_vgl & + (context, X, R, lmax, n, L, ldl, VGL, ldv) & + bind(C) result(info) + + use, intrinsic :: iso_c_binding + implicit none + + integer (c_int64_t) , intent(in) , value :: context + real (c_double ) , intent(in) :: X(3) + real (c_double ) , intent(in) :: R(3) + integer (c_int32_t) , intent(in) , value :: lmax + integer (c_int64_t) , intent(inout) :: n + integer (c_int32_t) , intent(out) :: L(ldl,n) + integer (c_int64_t) , intent(in) , value :: ldl + real (c_double ) , intent(out) :: VGL(ldv,5) + integer (c_int64_t) , intent(in) , value :: ldv + + integer(c_int32_t), external :: qmckl_ao_polynomial_transp_vgl_f + info = qmckl_ao_polynomial_transp_vgl_f & + (context, X, R, lmax, n, L, ldl, VGL, ldv) + + end function qmckl_ao_polynomial_transp_vgl + #+end_src + + #+CALL: generate_f_interface(table=qmckl_ao_polynomial_vgl_args,rettyp=get_value("FRetType"),fname=get_value("Name")) + + #+RESULTS: + #+begin_src f90 :tangle (eval fh_func) :comments org :exports none + interface + integer(c_int32_t) function qmckl_ao_polynomial_transp_vgl & + (context, X, R, lmax, n, L, ldl, VGL, ldv) & + bind(C) + use, intrinsic :: iso_c_binding + import + implicit none + + integer (c_int64_t) , intent(in) , value :: context + real (c_double ) , intent(in) :: X(3) + real (c_double ) , intent(in) :: R(3) + integer (c_int32_t) , intent(in) , value :: lmax + integer (c_int64_t) , intent(inout) :: n + integer (c_int32_t) , intent(out) :: L(ldl,n) + integer (c_int64_t) , intent(in) , value :: ldl + real (c_double ) , intent(out) :: VGL(ldv,5) + integer (c_int64_t) , intent(in) , value :: ldv + + end function qmckl_ao_polynomial_transp_vgl + end interface + #+end_src + *** Test :noexport: #+begin_src f90 :tangle (eval f_test) @@ -4348,11 +4551,12 @@ end function test_qmckl_ao_polynomial_vgl | ~nucleus_max_ang_mom~ | ~int32_t[nucl_num]~ | in | Maximum angular momentum per nucleus | | ~shell_ang_mom~ | ~int32_t[shell_num]~ | in | Angular momentum of each shell | | ~ao_factor~ | ~double[ao_num]~ | in | Normalization factor of the AOs | - | ~shell_vgl~ | ~double[5][point_num][shell_num]~ | in | Value, gradients and Laplacian of the shells | - | ~ao_vgl~ | ~double[5][point_num][ao_num]~ | out | Value, gradients and Laplacian of the AOs | + | ~shell_vgl~ | ~double[point_num][5][shell_num]~ | in | Value, gradients and Laplacian of the shells | + | ~ao_vgl~ | ~double[point_num][5][ao_num]~ | out | Value, gradients and Laplacian of the AOs | +** Unoptimized version #+begin_src f90 :comments org :tangle (eval f) :noweb yes -integer function qmckl_compute_ao_vgl_f(context, & +integer function qmckl_compute_ao_vgl_doc_f(context, & ao_num, shell_num, point_num, nucl_num, & coord, nucl_coord, nucleus_index, nucleus_shell_num, & nucleus_range, nucleus_max_ang_mom, shell_ang_mom, & @@ -4373,12 +4577,12 @@ integer function qmckl_compute_ao_vgl_f(context, & integer , intent(in) :: nucleus_max_ang_mom(nucl_num) integer , intent(in) :: shell_ang_mom(shell_num) double precision , intent(in) :: ao_factor(ao_num) - double precision , intent(in) :: shell_vgl(shell_num,point_num,5) - double precision , intent(out) :: ao_vgl(ao_num,point_num,5) + double precision , intent(in) :: shell_vgl(shell_num,5,point_num) + double precision , intent(out) :: ao_vgl(ao_num,5,point_num) double precision :: e_coord(3), n_coord(3) integer*8 :: n_poly - integer :: l, il, k, m, n + integer :: l, il, k integer*8 :: ipoint, inucl, ishell integer*8 :: ishell_start, ishell_end integer :: lstart(0:20) @@ -4388,102 +4592,229 @@ integer function qmckl_compute_ao_vgl_f(context, & double precision, allocatable :: poly_vgl(:,:) integer , allocatable :: powers(:,:) - integer , allocatable :: kil(:), knucl(:), kshell(:) - allocate(poly_vgl(8,ao_num), powers(8,ao_num)) - allocate(kil(ao_num), kshell(ao_num), knucl(nucl_num+1)) + allocate(poly_vgl(5,ao_num), powers(3,ao_num)) ! Pre-computed data - - k=1 - do inucl=1,nucl_num - knucl(inucl) = k - ishell_start = nucleus_index(inucl) + 1 - ishell_end = nucleus_index(inucl) + nucleus_shell_num(inucl) - do ishell = ishell_start, ishell_end - l = shell_ang_mom(ishell) - m = l*(l+1)*(l+2)/6 +1 - n = (l+1)*(l+2)*(l+3)/6 - do il = m, n - kil(k) = il - kshell(k) = ishell - k = k+1 - end do - end do + do l=0,20 + lstart(l) = l*(l+1)*(l+2)/6 +1 end do - knucl(nucl_num+1) = ao_num+1 - info = QMCKL_SUCCESS ! Don't compute polynomials when the radial part is zero. - ! TODO : Use numerical precision here cutoff = -dlog(1.d-12) do ipoint = 1, point_num e_coord(1) = coord(ipoint,1) e_coord(2) = coord(ipoint,2) e_coord(3) = coord(ipoint,3) - - ! Express the radial part in the AO basis - + k=1 do inucl=1,nucl_num - n_coord(1) = nucl_coord(inucl,1) - n_coord(2) = nucl_coord(inucl,2) - n_coord(3) = nucl_coord(inucl,3) + n_coord(1) = nucl_coord(inucl,1) + n_coord(2) = nucl_coord(inucl,2) + n_coord(3) = nucl_coord(inucl,3) ! Test if the point is in the range of the nucleus x = e_coord(1) - n_coord(1) y = e_coord(2) - n_coord(2) z = e_coord(3) - n_coord(3) - r2 = x*x + y*y + z*z + r2 = x*x + z*z + z*z if (r2 > cutoff*nucleus_range(inucl)) then - do k = knucl(inucl), knucl(inucl+1)-1 - ao_vgl(k,ipoint,1) = 0.d0 - ao_vgl(k,ipoint,2) = 0.d0 - ao_vgl(k,ipoint,3) = 0.d0 - ao_vgl(k,ipoint,4) = 0.d0 - ao_vgl(k,ipoint,5) = 0.d0 - end do cycle end if - ! Compute polynomials + ! Compute polynomials info = qmckl_ao_polynomial_vgl_f(context, e_coord, n_coord, & - nucleus_max_ang_mom(inucl), n_poly, powers, 8_8, poly_vgl, 8_8) + nucleus_max_ang_mom(inucl), n_poly, powers, 3_8, & + poly_vgl, 5_8) - do k = knucl(inucl), knucl(inucl+1)-1 - y = shell_vgl(kshell(k),ipoint,1) * ao_factor(k) + ! Loop over shells + ishell_start = nucleus_index(inucl) + 1 + ishell_end = nucleus_index(inucl) + nucleus_shell_num(inucl) + do ishell = ishell_start, ishell_end + l = shell_ang_mom(ishell) + do il = lstart(l), lstart(l+1)-1 + ! Value + ao_vgl(k,1,ipoint) = & + poly_vgl(1,il) * shell_vgl(ishell,1,ipoint) * ao_factor(k) - ao_vgl(k,ipoint,1) = y * poly_vgl(1,kil(k)) - ao_vgl(k,ipoint,2) = y * poly_vgl(2,kil(k)) - ao_vgl(k,ipoint,3) = y * poly_vgl(3,kil(k)) - ao_vgl(k,ipoint,4) = y * poly_vgl(4,kil(k)) - ao_vgl(k,ipoint,5) = y * poly_vgl(5,kil(k)) + ! Grad_x + ao_vgl(k,2,ipoint) = ( & + poly_vgl(2,il) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(1,il) * shell_vgl(ishell,2,ipoint) & + ) * ao_factor(k) - x = poly_vgl(1,kil(k)) * ao_factor(k) + ! Grad_y + ao_vgl(k,3,ipoint) = ( & + poly_vgl(3,il) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(1,il) * shell_vgl(ishell,3,ipoint) & + ) * ao_factor(k) - ao_vgl(k,ipoint,2) = ao_vgl(k,ipoint,2) + x * shell_vgl(kshell(k),ipoint,2) - ao_vgl(k,ipoint,3) = ao_vgl(k,ipoint,3) + x * shell_vgl(kshell(k),ipoint,3) - ao_vgl(k,ipoint,4) = ao_vgl(k,ipoint,4) + x * shell_vgl(kshell(k),ipoint,4) - ao_vgl(k,ipoint,5) = ao_vgl(k,ipoint,5) + x * shell_vgl(kshell(k),ipoint,5) + ! Grad_z + ao_vgl(k,4,ipoint) = ( & + poly_vgl(4,il) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(1,il) * shell_vgl(ishell,4,ipoint) & + ) * ao_factor(k) - ao_vgl(k,ipoint,5) = ao_vgl(k,ipoint,5) + & - (ao_factor(k) + ao_factor(k)) * (& - poly_vgl(2,kil(k)) * shell_vgl(kshell(k),ipoint,2) + & - poly_vgl(3,kil(k)) * shell_vgl(kshell(k),ipoint,3) + & - poly_vgl(4,kil(k)) * shell_vgl(kshell(k),ipoint,4) ) + ! Lapl_z + ao_vgl(k,5,ipoint) = ( & + poly_vgl(5,il) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(1,il) * shell_vgl(ishell,5,ipoint) + & + 2.d0 * ( & + poly_vgl(2,il) * shell_vgl(ishell,2,ipoint) + & + poly_vgl(3,il) * shell_vgl(ishell,3,ipoint) + & + poly_vgl(4,il) * shell_vgl(ishell,4,ipoint) ) & + ) * ao_factor(k) + + k = k+1 + end do end do - end do end do - deallocate(poly_vgl, powers, kshell, kil, knucl) -end function qmckl_compute_ao_vgl_f + deallocate(poly_vgl, powers) +end function qmckl_compute_ao_vgl_doc_f #+end_src +** HPC version + #+begin_src f90 :comments org :tangle (eval f) :noweb yes +integer function qmckl_compute_ao_vgl_hpc_f(context, & + ao_num, shell_num, point_num, nucl_num, & + coord, nucl_coord, nucleus_index, nucleus_shell_num, & + nucleus_range, nucleus_max_ang_mom, shell_ang_mom, & + ao_factor, shell_vgl, ao_vgl) & + result(info) + use qmckl + implicit none + integer(qmckl_context), intent(in) :: context + integer*8 , intent(in) :: ao_num + integer*8 , intent(in) :: shell_num + integer*8 , intent(in) :: point_num + integer*8 , intent(in) :: nucl_num + double precision , intent(in) :: coord(point_num,3) + double precision , intent(in) :: nucl_coord(nucl_num,3) + integer*8 , intent(in) :: nucleus_index(nucl_num) + integer*8 , intent(in) :: nucleus_shell_num(nucl_num) + double precision , intent(in) :: nucleus_range(nucl_num) + integer , intent(in) :: nucleus_max_ang_mom(nucl_num) + integer , intent(in) :: shell_ang_mom(shell_num) + double precision , intent(in) :: ao_factor(ao_num) + double precision , intent(in) :: shell_vgl(shell_num,5,point_num) + double precision , intent(out) :: ao_vgl(ao_num,5,point_num) + + double precision :: e_coord(3), n_coord(3) + integer*8 :: n_poly + integer :: l, il, k + integer*8 :: ipoint, inucl, ishell + integer*8 :: ishell_start, ishell_end + integer :: lstart(0:20) + double precision :: x, y, z, r2 + double precision :: cutoff + integer, external :: qmckl_ao_polynomial_transp_vgl_f + + double precision, allocatable :: poly_vgl(:,:) + integer , allocatable :: powers(:,:) + + allocate(poly_vgl(ao_num,5), powers(3,ao_num)) + + ! Pre-computed data + do l=0,20 + lstart(l) = l*(l+1)*(l+2)/6 +1 + end do + + info = QMCKL_SUCCESS + + ! Don't compute polynomials when the radial part is zero. + cutoff = -dlog(1.d-12) + + do ipoint = 1, point_num + e_coord(1) = coord(ipoint,1) + e_coord(2) = coord(ipoint,2) + e_coord(3) = coord(ipoint,3) + k=1 + do inucl=1,nucl_num + n_coord(1) = nucl_coord(inucl,1) + n_coord(2) = nucl_coord(inucl,2) + n_coord(3) = nucl_coord(inucl,3) + + ! Test if the point is in the range of the nucleus + x = e_coord(1) - n_coord(1) + y = e_coord(2) - n_coord(2) + z = e_coord(3) - n_coord(3) + + r2 = x*x + z*z + z*z + + if (r2 > cutoff*nucleus_range(inucl)) then + cycle + end if + + ! Compute polynomials + info = qmckl_ao_polynomial_transp_vgl_f(context, e_coord, n_coord, & + nucleus_max_ang_mom(inucl), n_poly, powers, 3_8, & + poly_vgl, int(ao_num,8)) + + ! Loop over shells + ishell_start = nucleus_index(inucl) + 1 + ishell_end = nucleus_index(inucl) + nucleus_shell_num(inucl) + do ishell = ishell_start, ishell_end + l = shell_ang_mom(ishell) + if (shell_vgl(ishell,1,ipoint) /= 0.d0) then + do il = lstart(l), lstart(l+1)-1 + ! Value + ao_vgl(k,1,ipoint) = & + poly_vgl(il,1) * shell_vgl(ishell,1,ipoint) * ao_factor(k) + + ! Grad_x + ao_vgl(k,2,ipoint) = ( & + poly_vgl(il,2) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(il,1) * shell_vgl(ishell,2,ipoint) & + ) * ao_factor(k) + + ! Grad_y + ao_vgl(k,3,ipoint) = ( & + poly_vgl(il,3) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(il,1) * shell_vgl(ishell,3,ipoint) & + ) * ao_factor(k) + + ! Grad_z + ao_vgl(k,4,ipoint) = ( & + poly_vgl(il,4) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(il,1) * shell_vgl(ishell,4,ipoint) & + ) * ao_factor(k) + + ! Lapl_z + ao_vgl(k,5,ipoint) = ( & + poly_vgl(il,5) * shell_vgl(ishell,1,ipoint) + & + poly_vgl(il,1) * shell_vgl(ishell,5,ipoint) + & + 2.d0 * ( & + poly_vgl(il,2) * shell_vgl(ishell,2,ipoint) + & + poly_vgl(il,3) * shell_vgl(ishell,3,ipoint) + & + poly_vgl(il,4) * shell_vgl(ishell,4,ipoint) ) & + ) * ao_factor(k) + k = k+1 + end do + else + do il = lstart(l), lstart(l+1)-1 + ao_vgl(k,1,ipoint) = 0.d0 + ao_vgl(k,2,ipoint) = 0.d0 + ao_vgl(k,3,ipoint) = 0.d0 + ao_vgl(k,4,ipoint) = 0.d0 + ao_vgl(k,5,ipoint) = 0.d0 + k = k+1 + end do + end if + end do + end do + end do + + deallocate(poly_vgl, powers) +end function qmckl_compute_ao_vgl_hpc_f + #+end_src + +** Interfaces # #+CALL: generate_c_header(table=qmckl_ao_vgl_args,rettyp=get_value("CRetType"),fname="qmckl_compute_ao_vgl")) # (Commented because the header needs to go into h_private_func @@ -4545,11 +4876,16 @@ end function qmckl_compute_ao_vgl_f integer (c_int32_t) , intent(in) :: nucleus_max_ang_mom(nucl_num) integer (c_int32_t) , intent(in) :: shell_ang_mom(shell_num) real (c_double ) , intent(in) :: ao_factor(ao_num) - real (c_double ) , intent(in) :: shell_vgl(shell_num,point_num,5) - real (c_double ) , intent(out) :: ao_vgl(ao_num,point_num,5) + real (c_double ) , intent(in) :: shell_vgl(shell_num,5,point_num) + real (c_double ) , intent(out) :: ao_vgl(ao_num,5,point_num) - integer(c_int32_t), external :: qmckl_compute_ao_vgl_f - info = qmckl_compute_ao_vgl_f & +#ifdef HAVE_HPC + integer(c_int32_t), external :: qmckl_compute_ao_vgl_hpc_f + info = qmckl_compute_ao_vgl_hpc_f & +#else + integer(c_int32_t), external :: qmckl_compute_ao_vgl_doc_f + info = qmckl_compute_ao_vgl_doc_f & +#endif (context, & ao_num, & shell_num, & @@ -4750,39 +5086,39 @@ rc = qmckl_set_electron_coord (context, 'N', elec_coord, walk_num*elec_num*3); assert(rc == QMCKL_SUCCESS); -double ao_vgl[5][elec_num][ao_num]; +double ao_vgl[elec_num][5][ao_num]; rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]), (int64_t) 5*elec_num*ao_num); assert (rc == QMCKL_SUCCESS); printf("\n"); -printf(" ao_vgl ao_vgl[0][26][219] %25.15e\n", ao_vgl[0][26][219]); -printf(" ao_vgl ao_vgl[1][26][219] %25.15e\n", ao_vgl[1][26][219]); -printf(" ao_vgl ao_vgl[0][26][220] %25.15e\n", ao_vgl[0][26][220]); -printf(" ao_vgl ao_vgl[1][26][220] %25.15e\n", ao_vgl[1][26][220]); -printf(" ao_vgl ao_vgl[0][26][221] %25.15e\n", ao_vgl[0][26][221]); -printf(" ao_vgl ao_vgl[1][26][221] %25.15e\n", ao_vgl[1][26][221]); -printf(" ao_vgl ao_vgl[0][26][222] %25.15e\n", ao_vgl[0][26][222]); -printf(" ao_vgl ao_vgl[1][26][222] %25.15e\n", ao_vgl[1][26][222]); -printf(" ao_vgl ao_vgl[0][26][223] %25.15e\n", ao_vgl[0][26][223]); -printf(" ao_vgl ao_vgl[1][26][223] %25.15e\n", ao_vgl[1][26][223]); -printf(" ao_vgl ao_vgl[0][26][224] %25.15e\n", ao_vgl[0][26][224]); -printf(" ao_vgl ao_vgl[1][26][224] %25.15e\n", ao_vgl[1][26][224]); +printf(" ao_vgl ao_vgl[0][26][219] %25.15e\n", ao_vgl[26][0][219]); +printf(" ao_vgl ao_vgl[1][26][219] %25.15e\n", ao_vgl[26][1][219]); +printf(" ao_vgl ao_vgl[0][26][220] %25.15e\n", ao_vgl[26][0][220]); +printf(" ao_vgl ao_vgl[1][26][220] %25.15e\n", ao_vgl[26][1][220]); +printf(" ao_vgl ao_vgl[0][26][221] %25.15e\n", ao_vgl[26][0][221]); +printf(" ao_vgl ao_vgl[1][26][221] %25.15e\n", ao_vgl[26][1][221]); +printf(" ao_vgl ao_vgl[0][26][222] %25.15e\n", ao_vgl[26][0][222]); +printf(" ao_vgl ao_vgl[1][26][222] %25.15e\n", ao_vgl[26][1][222]); +printf(" ao_vgl ao_vgl[0][26][223] %25.15e\n", ao_vgl[26][0][223]); +printf(" ao_vgl ao_vgl[1][26][223] %25.15e\n", ao_vgl[26][1][223]); +printf(" ao_vgl ao_vgl[0][26][224] %25.15e\n", ao_vgl[26][0][224]); +printf(" ao_vgl ao_vgl[1][26][224] %25.15e\n", ao_vgl[26][1][224]); printf("\n"); -assert( fabs(ao_vgl[0][26][219] - ( 1.020298798341620e-08)) < 1.e-14 ); -assert( fabs(ao_vgl[1][26][219] - (-4.928035238010602e-08)) < 1.e-14 ); -assert( fabs(ao_vgl[0][26][220] - ( 1.516643537739178e-08)) < 1.e-14 ); -assert( fabs(ao_vgl[1][26][220] - (-7.725221462603871e-08)) < 1.e-14 ); -assert( fabs(ao_vgl[0][26][221] - (-4.686370882518819e-09)) < 1.e-14 ); -assert( fabs(ao_vgl[1][26][221] - ( 2.387064067626827e-08)) < 1.e-14 ); -assert( fabs(ao_vgl[0][26][222] - ( 7.514816980753531e-09)) < 1.e-14 ); -assert( fabs(ao_vgl[1][26][222] - (-4.025889138635182e-08)) < 1.e-14 ); -assert( fabs(ao_vgl[0][26][223] - (-4.021908374204471e-09)) < 1.e-14 ); -assert( fabs(ao_vgl[1][26][223] - ( 2.154644255710413e-08)) < 1.e-14 ); -assert( fabs(ao_vgl[0][26][224] - ( 7.175045873560788e-10)) < 1.e-14 ); -assert( fabs(ao_vgl[1][26][224] - (-3.843864637762753e-09)) < 1.e-14 ); +assert( fabs(ao_vgl[26][0][219] - ( 1.020298798341620e-08)) < 1.e-14 ); +assert( fabs(ao_vgl[26][1][219] - (-4.928035238010602e-08)) < 1.e-14 ); +assert( fabs(ao_vgl[26][0][220] - ( 1.516643537739178e-08)) < 1.e-14 ); +assert( fabs(ao_vgl[26][1][220] - (-7.725221462603871e-08)) < 1.e-14 ); +assert( fabs(ao_vgl[26][0][221] - (-4.686370882518819e-09)) < 1.e-14 ); +assert( fabs(ao_vgl[26][1][221] - ( 2.387064067626827e-08)) < 1.e-14 ); +assert( fabs(ao_vgl[26][0][222] - ( 7.514816980753531e-09)) < 1.e-14 ); +assert( fabs(ao_vgl[26][1][222] - (-4.025889138635182e-08)) < 1.e-14 ); +assert( fabs(ao_vgl[26][0][223] - (-4.021908374204471e-09)) < 1.e-14 ); +assert( fabs(ao_vgl[26][1][223] - ( 2.154644255710413e-08)) < 1.e-14 ); +assert( fabs(ao_vgl[26][0][224] - ( 7.175045873560788e-10)) < 1.e-14 ); +assert( fabs(ao_vgl[26][1][224] - (-3.843864637762753e-09)) < 1.e-14 ); } #+end_src diff --git a/org/qmckl_mo.org b/org/qmckl_mo.org index 462828d..709fb31 100644 --- a/org/qmckl_mo.org +++ b/org/qmckl_mo.org @@ -90,11 +90,11 @@ int main() { Computed data: - |---------------+-------------------------+----------------------------------------------------------------------------------------| - |---------------+-------------------------+----------------------------------------------------------------------------------------| - | ~mo_vgl~ | ~[5][elec_num][mo_num]~ | Value, gradients, Laplacian of the MOs at electron positions | - | ~mo_vgl_date~ | ~uint64_t~ | Late modification date of Value, gradients, Laplacian of the MOs at electron positions | - |---------------+-------------------------+----------------------------------------------------------------------------------------| + |---------------+--------------------------+-------------------------------------------------------------------------------------| + |---------------+--------------------------+-------------------------------------------------------------------------------------| + | ~mo_vgl~ | ~[point_num][5][mo_num]~ | Value, gradients, Laplacian of the MOs at point positions | + | ~mo_vgl_date~ | ~uint64_t~ | Late modification date of Value, gradients, Laplacian of the MOs at point positions | + |---------------+--------------------------+-------------------------------------------------------------------------------------| ** Data structure @@ -388,7 +388,7 @@ qmckl_exit_code qmckl_get_mo_basis_vgl(qmckl_context context, double* const mo_v qmckl_context_struct* const ctx = (qmckl_context_struct* const) context; assert (ctx != NULL); - size_t sze = 5 * ctx->electron.num * ctx->mo_basis.mo_num; + size_t sze = 5 * ctx->point.num * ctx->mo_basis.mo_num; memcpy(mo_vgl, ctx->mo_basis.mo_vgl, sze * sizeof(double)); return QMCKL_SUCCESS; @@ -442,13 +442,6 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) NULL); } - if(!(ctx->electron.provided)) { - return qmckl_failwith( context, - QMCKL_NOT_PROVIDED, - "qmckl_electron", - NULL); - } - if (!ctx->mo_basis.provided) { return qmckl_failwith( context, QMCKL_NOT_PROVIDED, @@ -457,13 +450,13 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) } /* Compute if necessary */ - if (ctx->electron.coord_new_date > ctx->mo_basis.mo_vgl_date) { + if (ctx->point.date > ctx->mo_basis.mo_vgl_date) { /* Allocate array */ if (ctx->mo_basis.mo_vgl == NULL) { qmckl_memory_info_struct mem_info = qmckl_memory_info_struct_zero; - mem_info.size = 5 * ctx->electron.num * ctx->mo_basis.mo_num * sizeof(double); + mem_info.size = 5 * ctx->point.num * ctx->mo_basis.mo_num * sizeof(double); double* mo_vgl = (double*) qmckl_malloc(context, mem_info); if (mo_vgl == NULL) { @@ -478,7 +471,7 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) rc = qmckl_compute_mo_basis_vgl(context, ctx->ao_basis.ao_num, ctx->mo_basis.mo_num, - ctx->electron.num, + ctx->point.num, ctx->mo_basis.coefficient, ctx->ao_basis.ao_vgl, ctx->mo_basis.mo_vgl); @@ -504,85 +497,46 @@ qmckl_exit_code qmckl_provide_mo_vgl(qmckl_context context) | ~qmckl_context~ | ~context~ | in | Global state | | ~int64_t~ | ~ao_num~ | in | Number of AOs | | ~int64_t~ | ~mo_num~ | in | Number of MOs | - | ~int64_t~ | ~elec_num~ | in | Number of electrons | + | ~int64_t~ | ~point_num~ | in | Number of points | | ~double~ | ~coef_normalized[mo_num][ao_num]~ | in | AO to MO transformation matrix | - | ~double~ | ~ao_vgl[5][elec_num][ao_num]~ | in | Value, gradients and Laplacian of the AOs | - | ~double~ | ~mo_vgl[5][elec_num][mo_num]~ | out | Value, gradients and Laplacian of the MOs | + | ~double~ | ~ao_vgl[point_num][5][ao_num]~ | in | Value, gradients and Laplacian of the AOs | + | ~double~ | ~mo_vgl[point_num][5][mo_num]~ | out | Value, gradients and Laplacian of the MOs | #+begin_src f90 :comments org :tangle (eval f) :noweb yes integer function qmckl_compute_mo_basis_vgl_f(context, & - ao_num, mo_num, elec_num, & + ao_num, mo_num, point_num, & coef_normalized, ao_vgl, mo_vgl) & result(info) use qmckl implicit none integer(qmckl_context), intent(in) :: context integer*8 , intent(in) :: ao_num, mo_num - integer*8 , intent(in) :: elec_num - double precision , intent(in) :: ao_vgl(ao_num,elec_num,5) + integer*8 , intent(in) :: point_num + double precision , intent(in) :: ao_vgl(ao_num,5,point_num) double precision , intent(in) :: coef_normalized(ao_num,mo_num) - double precision , intent(out) :: mo_vgl(mo_num,elec_num,5) + double precision , intent(out) :: mo_vgl(mo_num,5,point_num) character :: TransA, TransB - double precision,dimension(:,:),allocatable :: mo_vgl_big - double precision,dimension(:,:),allocatable :: ao_vgl_big - !double precision,dimension(:,:),allocatable :: coef_trans - !double precision,dimension(:),allocatable :: coef_all double precision :: alpha, beta - integer :: info_qmckl_dgemm_value - integer*8 :: M, N, K, LDA, LDB, LDC, i,j, idx - - integer*8 :: inucl, iprim, iwalk, ielec, ishell - double precision :: x, y, z, two_a, ar2, r2, v, cutoff - - allocate(mo_vgl_big(mo_num,elec_num*5)) - allocate(ao_vgl_big(ao_num,elec_num*5)) - !allocate(coef_all(mo_num*ao_num)) - !allocate(coef_trans(mo_num,ao_num)) + integer*8 :: M, N, K, LDA, LDB, LDC, i,j TransA = 'T' TransB = 'N' - alpha = 1.0d0 - beta = 0.0d0 + M = mo_num + N = point_num*5_8 + K = int(ao_num,8) + alpha = 1.d0 + beta = 0.d0 + LDA = size(coef_normalized,1) + LDB = size(ao_vgl,1) + LDC = size(mo_vgl,1) + + info = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, & + coef_normalized, int(size(coef_normalized,1),8), & + ao_vgl, int(size(ao_vgl,1),8), beta, & + mo_vgl,LDC) info = QMCKL_SUCCESS - info_qmckl_dgemm_value = QMCKL_SUCCESS - - ! Don't compute exponentials when the result will be almost zero. - ! TODO : Use numerical precision here - cutoff = -dlog(1.d-15) - M = mo_num - N = elec_num*5 - K = ao_num * 1_8 - LDA = size(coef_normalized,1) - idx = 0 - !do j = 1,ao_num - !do i = 1,mo_num - ! idx = idx + 1 - ! coef_all(idx) = coef_normalized(i,j) - !end do - !end do - !idx = 0 - !do j = 1,mo_num - !do i = 1,ao_num - ! idx = idx + 1 - ! coef_trans(j,i) = coef_all(idx) - !end do - !end do - - ao_vgl_big = reshape(ao_vgl(:, :, :),(/ao_num, elec_num*5_8/)) - LDB = size(ao_vgl_big,1) - LDC = size(mo_vgl_big,1) - - info = qmckl_dgemm(context,TransA, TransB, M, N, K, alpha, & - coef_normalized,size(coef_normalized,1)*1_8, & - ao_vgl_big, size(ao_vgl_big,1)*1_8, & - beta, & - mo_vgl_big,LDC) - mo_vgl = reshape(mo_vgl_big,(/mo_num,elec_num,5_8/)) - - deallocate(mo_vgl_big) - deallocate(ao_vgl_big) end function qmckl_compute_mo_basis_vgl_f #+end_src @@ -595,7 +549,7 @@ end function qmckl_compute_mo_basis_vgl_f const qmckl_context context, const int64_t ao_num, const int64_t mo_num, - const int64_t elec_num, + const int64_t point_num, const double* coef_normalized, const double* ao_vgl, double* const mo_vgl ); @@ -607,7 +561,7 @@ end function qmckl_compute_mo_basis_vgl_f #+RESULTS: #+begin_src f90 :tangle (eval f) :comments org :exports none integer(c_int32_t) function qmckl_compute_mo_basis_vgl & - (context, ao_num, mo_num, elec_num, coef_normalized, ao_vgl, mo_vgl) & + (context, ao_num, mo_num, point_num, coef_normalized, ao_vgl, mo_vgl) & bind(C) result(info) use, intrinsic :: iso_c_binding @@ -616,14 +570,14 @@ end function qmckl_compute_mo_basis_vgl_f integer (c_int64_t) , intent(in) , value :: context integer (c_int64_t) , intent(in) , value :: ao_num integer (c_int64_t) , intent(in) , value :: mo_num - integer (c_int64_t) , intent(in) , value :: elec_num + integer (c_int64_t) , intent(in) , value :: point_num real (c_double ) , intent(in) :: coef_normalized(ao_num,mo_num) - real (c_double ) , intent(in) :: ao_vgl(ao_num,elec_num,5) - real (c_double ) , intent(out) :: mo_vgl(mo_num,elec_num,5) + real (c_double ) , intent(in) :: ao_vgl(ao_num,5,point_num) + real (c_double ) , intent(out) :: mo_vgl(mo_num,5,point_num) integer(c_int32_t), external :: qmckl_compute_mo_basis_vgl_f info = qmckl_compute_mo_basis_vgl_f & - (context, ao_num, mo_num, elec_num, coef_normalized, ao_vgl, mo_vgl) + (context, ao_num, mo_num, point_num, coef_normalized, ao_vgl, mo_vgl) end function qmckl_compute_mo_basis_vgl #+end_src @@ -799,9 +753,9 @@ assert(rc == QMCKL_SUCCESS); assert(qmckl_ao_basis_provided(context)); -double ao_vgl[5][walk_num][elec_num][chbrclf_ao_num]; +double ao_vgl[walk_num*elec_num][5][chbrclf_ao_num]; -rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0][0]), +rc = qmckl_get_ao_basis_ao_vgl(context, &(ao_vgl[0][0][0]), (int64_t) 5*walk_num*elec_num*chbrclf_ao_num); assert (rc == QMCKL_SUCCESS); @@ -817,7 +771,7 @@ assert (rc == QMCKL_SUCCESS); assert(qmckl_mo_basis_provided(context)); -double mo_vgl[5][elec_num][chbrclf_mo_num]; +double mo_vgl[walk_num*elec_num][5][chbrclf_mo_num]; rc = qmckl_get_mo_basis_vgl(context, &(mo_vgl[0][0][0])); assert (rc == QMCKL_SUCCESS); @@ -863,18 +817,18 @@ assert (rc == QMCKL_SUCCESS); printf("\n"); -printf(" mo_vgl mo_vgl[0][26][219] %25.15e\n", mo_vgl[0][2][3]); -printf(" mo_vgl mo_vgl[1][26][219] %25.15e\n", mo_vgl[1][2][3]); -printf(" mo_vgl mo_vgl[0][26][220] %25.15e\n", mo_vgl[0][2][3]); -printf(" mo_vgl mo_vgl[1][26][220] %25.15e\n", mo_vgl[1][2][3]); -printf(" mo_vgl mo_vgl[0][26][221] %25.15e\n", mo_vgl[0][2][3]); -printf(" mo_vgl mo_vgl[1][26][221] %25.15e\n", mo_vgl[1][2][3]); -printf(" mo_vgl mo_vgl[0][26][222] %25.15e\n", mo_vgl[0][2][3]); -printf(" mo_vgl mo_vgl[1][26][222] %25.15e\n", mo_vgl[1][2][3]); -printf(" mo_vgl mo_vgl[0][26][223] %25.15e\n", mo_vgl[0][2][3]); -printf(" mo_vgl mo_vgl[1][26][223] %25.15e\n", mo_vgl[1][2][3]); -printf(" mo_vgl mo_vgl[0][26][224] %25.15e\n", mo_vgl[0][2][3]); -printf(" mo_vgl mo_vgl[1][26][224] %25.15e\n", mo_vgl[1][2][3]); +printf(" mo_vgl mo_vgl[0][26][219] %25.15e\n", mo_vgl[2][0][3]); +printf(" mo_vgl mo_vgl[1][26][219] %25.15e\n", mo_vgl[2][1][3]); +printf(" mo_vgl mo_vgl[0][26][220] %25.15e\n", mo_vgl[2][0][3]); +printf(" mo_vgl mo_vgl[1][26][220] %25.15e\n", mo_vgl[2][1][3]); +printf(" mo_vgl mo_vgl[0][26][221] %25.15e\n", mo_vgl[2][0][3]); +printf(" mo_vgl mo_vgl[1][26][221] %25.15e\n", mo_vgl[2][1][3]); +printf(" mo_vgl mo_vgl[0][26][222] %25.15e\n", mo_vgl[2][0][3]); +printf(" mo_vgl mo_vgl[1][26][222] %25.15e\n", mo_vgl[2][1][3]); +printf(" mo_vgl mo_vgl[0][26][223] %25.15e\n", mo_vgl[2][0][3]); +printf(" mo_vgl mo_vgl[1][26][223] %25.15e\n", mo_vgl[2][1][3]); +printf(" mo_vgl mo_vgl[0][26][224] %25.15e\n", mo_vgl[2][0][3]); +printf(" mo_vgl mo_vgl[1][26][224] %25.15e\n", mo_vgl[2][1][3]); printf("\n"); } diff --git a/tools/Building.org b/tools/Building.org deleted file mode 100644 index 2cfa91c..0000000 --- a/tools/Building.org +++ /dev/null @@ -1,641 +0,0 @@ -#+TITLE: Building tools -#+STARTUP: indent overview -#+PROPERTY: header-args: :comments both - -This file contains all the tools needed to build the QMCkl library. - -* Helper functions - #+NAME: header - #+begin_src sh :tangle no :exports none :output none -echo "This file was created by tools/Building.org" - #+end_src - - #+NAME: check-src - #+begin_src bash -if [[ $(basename ${PWD}) != "src" ]] ; then - echo "This script needs to be run in the src directory" - exit -1 -fi - #+end_src - - #+NAME: url-issues - : https://github.com/trex-coe/qmckl/issues - - #+NAME: url-web - : https://trex-coe.github.io/qmckl - - #+NAME: license - #+begin_example -BSD 3-Clause License - -Copyright (c) 2020, TREX Center of Excellence -All rights reserved. - -Redistribution and use in source and binary forms, with or without -modification, are permitted provided that the following conditions are met: - -1. Redistributions of source code must retain the above copyright notice, this - list of conditions and the following disclaimer. - -2. Redistributions in binary form must reproduce the above copyright notice, - this list of conditions and the following disclaimer in the documentation - and/or other materials provided with the distribution. - -3. Neither the name of the copyright holder nor the names of its - contributors may be used to endorse or promote products derived from - this software without specific prior written permission. - -THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" -AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE -IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE -DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE -FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL -DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR -SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER -CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, -OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE -OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. - - #+end_example - -* Makefiles -** Makefile.in -:PROPERTIES: -:header-args: :tangle ../src/Makefile.in :noweb yes :comments org -:END: - -This is the main Makefile invoked by the ~make~ command at the root -of the package. To compile the sources, it calls the =Makefile= -located in the =src= directory. This Makefile creates the source -file from the org-mode file, as well as a Makefile, -=Makefile.generated=, dedicated to the compilation of the sources. - -*** Header - -We want the Makefile to be POSIX-compliant, such that it works not -only with GNU Make. - -#+begin_src makefile -# <> - -.POSIX: -#+end_src - -*** Compiler options - -Compiler variables are obtained from the configure script (see =configure.ac=) - -#+begin_src makefile -CC = @CC@ -FC = @FC@ -CFLAGS = @CFLAGS@ -FCFLAGS = @FCFLAGS@ -LDFLAGS = @LDFLAGS@ -DEFS = @DEFS@ - -#+end_src - -*** Variables - -#+begin_src makefile -HAS_CPPCHECK = @HAS_CPPCHECK@ - -# VPATH-related substitution variables -srcdir = @srcdir@ -VPATH = @srcdir@ - -top_srcdir=$(srcdir)/.. -shared_lib=$(top_srcdir)/lib/libqmckl.so -static_lib=$(top_srcdir)/lib/libqmckl.a -qmckl_h=$(top_srcdir)/include/qmckl.h -qmckl_f=$(top_srcdir)/share/qmckl/fortran/qmckl_f.f90 - -export CC CFLAGS DEFS FC FCFLAGS LIBS top_srcdir - -ORG_SOURCE_FILES=$(wildcard $(srcdir)/*.org) -C_SOURCE_FILES=$(patsubst %.org,%.c,$(ORG_SOURCE_FILES)) -INCLUDE=-I$(top_srcdir)/include/ -#+end_src - -*** Rules - -The source files are created during the generation of the file ~Makefile.generated~. -The Makefile.generated is the one that will be distributed with the library. - -#+begin_src makefile -.PHONY: clean shared static doc all check install uninstall -.SECONDARY: # Needed to keep the produced C and Fortran files - -$(shared_lib) $(static_lib): $(qmckl_h) $(qmckl_f) Makefile.generated - $(MAKE) -f Makefile.generated $@ - -install uninstall: Makefile.generated - $(MAKE) -f Makefile.generated $@ - -$(qmckl_f) $(qmckl_h): Makefile.generated - $(top_srcdir)/tools/build_qmckl_h.sh - -shared: $(shared_lib) -static: $(static_lib) -all: shared static doc check - -check: $(static_lib) - $(MAKE) -f Makefile.generated check - -ifeq ($(HAS_CPPCHECK),1) -cppcheck: - cppcheck \ - --addon=cert \ - --enable=warning,style,performance,portability,information \ - qmckl_*.c -endif - -doc: $(ORG_SOURCE_FILES) - $(top_srcdir)/tools/build_doc.sh - -clean: - - $(MAKE) -f Makefile.generated clean - - $(RM) test_qmckl_* test_qmckl.c \ - $(qmckl_h) $(qmckl_f) \ - qmckl_*.f90 qmckl_*.c qmckl_*.h \ - Makefile.generated *.html *.txt - -veryclean: clean FORCE - - $(RM) $(top_srcdir)/share/doc/qmckl/html/*.html \ - $(top_srcdir)/share/doc/qmckl/text/*.txt - -Makefile.generated.in: Makefile $(top_srcdir)/tools/create_makefile.sh $(ORG_SOURCE_FILES) $(top_srcdir)/tools/Building.org - $(top_srcdir)/tools/create_makefile.sh - -Makefile.generated: Makefile.generated.in - cd .. ; ./config.status - -.SUFFIXES: .org .c - -.org.c: - $(top_srcdir)/tools/tangle.sh $< - -#+end_src - -** Script to generate auto-generated Makefile - :PROPERTIES: - :header-args: :tangle create_makefile.sh :noweb yes :shebang #!/bin/bash :comments org - :END: - - This script generates the Makefile that compiles the library. - The ~OUTPUT~ variable contains the name of the generated Makefile,typically - =Makefile.generated=. - - #+begin_src bash -# <> - -<> - -OUTPUT=Makefile.generated.in - #+end_src - - We start by tangling all the org-mode files. - - #+begin_src bash -${top_srcdir}/tools/tangle.sh *.org -${top_srcdir}/tools/build_qmckl_h.sh - #+end_src - - Then we create the list of ~*.o~ files to be created, for library - functions: - - #+begin_src bash -OBJECTS="qmckl_f.o" -for i in $(ls qmckl_*.c qmckl_*f.f90) ; do - FILE=${i%.*} - OBJECTS+=" ${FILE}.o" -done >> $OUTPUT - #+end_src - - for tests in C: - - #+begin_src bash -TESTS="" -for i in $(ls test_qmckl_*.c) ; do - FILE=${i%.c} - TESTS+=" ${FILE}.o" -done >> $OUTPUT - #+end_src - - and for tests in Fortran: - - #+begin_src bash -TESTS_F="" -for i in $(ls test_qmckl_*_f.f90) ; do - FILE=${i%.f90} - TESTS_F+=" ${FILE}.o" -done >> $OUTPUT - #+end_src - - Finally, we append the variables to the Makefile - - #+begin_src bash :noweb yes -cat << EOF > ${OUTPUT} -.POSIX: -.SUFFIXES: - -package = @PACKAGE_TARNAME@ -version = @PACKAGE_VERSION@ - -# VPATH-related substitution variables -srcdir = @srcdir@ -VPATH = @srcdir@ - -prefix = @prefix@ - -CC = @CC@ -DEFS = @DEFS@ -CFLAGS = @CFLAGS@ -I\$(top_srcdir)/munit/ -I\$(top_srcdir)/include -I. -CPPFLAGS = @CPPFLAGS@ -LIBS = @LIBS@ - -FC = @FC@ -FCFLAGS= @FCFLAGS@ - -OBJECT_FILES=$OBJECTS - -TESTS = $TESTS -TESTS_F = $TESTS_F - -LIBS = @LIBS@ -FCLIBS = @FCLIBS@ -EOF - -export -echo ' -<> -' >> ${OUTPUT} - - #+end_src - -and the rules: - -#+NAME: rules - #+begin_src makefile :tangle no -top_srcdir=$(srcdir)/.. -shared_lib=$(top_srcdir)/lib/libqmckl.so -static_lib=$(top_srcdir)/lib/libqmckl.a -qmckl_h=$(top_srcdir)/include/qmckl.h -qmckl_f=$(top_srcdir)/share/qmckl/fortran/qmckl_f.f90 -munit=$(top_srcdir)/munit/munit.c - -datarootdir=$(prefix)/share -datadir=$(datarootdir) -docdir=$(datarootdir)/doc/$(package) -htmldir=$(docdir)/html -libdir=$(prefix)/lib -includedir=$(prefix)/include -fortrandir=$(datarootdir)/$(package)/fortran - - -shared: $(shared_lib) -static: $(static_lib) - - -all: shared static - -$(shared_lib): $(OBJECT_FILES) - $(CC) -shared $(OBJECT_FILES) -o $(shared_lib) - -$(static_lib): $(OBJECT_FILES) - $(AR) rcs $(static_lib) $(OBJECT_FILES) - - -# Test - -qmckl_f.o: $(qmckl_f) - $(FC) $(FCFLAGS) -c $(qmckl_f) -o $@ - -test_qmckl: test_qmckl.c $(qmckl_h) $(static_lib) $(TESTS) $(TESTS_F) - $(CC) $(CFLAGS) $(CPPFLAGS) $(DEFS) $(munit) $(TESTS) $(TESTS_F) \ - $(static_lib) $(LIBS) $(FCLIBS) test_qmckl.c -o $@ - -test_qmckl_shared: test_qmckl.c $(qmckl_h) $(shared_lib) $(TESTS) $(TESTS_F) - $(CC) $(CFLAGS) $(CPPFLAGS) $(DEFS) \ - -Wl,-rpath,$(top_srcdir)/lib -L$(top_srcdir)/lib $(munit) $(TESTS) \ - $(TESTS_F) -lqmckl $(LIBS) $(FCLIBS) test_qmckl.c -o $@ - -check: test_qmckl test_qmckl_shared - ./test_qmckl - -clean: - $(RM) -- *.o *.mod $(shared_lib) $(static_lib) test_qmckl - - - - -install: - install -d $(DESTDIR)$(prefix)/lib - install -d $(DESTDIR)$(prefix)/include - install -d $(DESTDIR)$(prefix)/share/qmckl/fortran - install -d $(DESTDIR)$(prefix)/share/doc/qmckl/html/ - install -d $(DESTDIR)$(prefix)/share/doc/qmckl/text/ - install $(shared_lib) $(DESTDIR)$(libdir)/ - install $(static_lib) $(DESTDIR)$(libdir)/ - install $(qmckl_h) $(DESTDIR)$(includedir) - install $(qmckl_f) $(DESTDIR)$(fortrandir) - install $(top_srcdir)/share/doc/qmckl/html/*.html $(DESTDIR)$(docdir)/html/ - install $(top_srcdir)/share/doc/qmckl/html/*.css $(DESTDIR)$(docdir)/html/ - install $(top_srcdir)/share/doc/qmckl/text/*.txt $(DESTDIR)$(docdir)/text/ - -uninstall: - rm $(DESTDIR)$(libdir)/libqmckl.so - rm $(DESTDIR)$(libdir)/libqmckl.a - rm $(DESTDIR)$(includedir)/qmckl.h - rm -rf $(DESTDIR)$(datarootdir)/$(package) - rm -rf $(DESTDIR)$(docdir) - -.SUFFIXES: .c .f90 .o - -.c.o: - $(CC) $(CFLAGS) $(CPPFLAGS) $(DEFS) -c $*.c -o $*.o - -.f90.o: qmckl_f.o - $(FC) $(FCFLAGS) -c $*.f90 -o $*.o - -.PHONY: check cppcheck clean all - #+end_src - -* Script to tangle the org-mode files - :PROPERTIES: - :header-args: :tangle tangle.sh :noweb yes :shebang #!/bin/bash :comments org - :END: - - #+begin_src bash -# <> - -<> - #+end_src - - This file needs to be run from the QMCKL =src= directory. - - It tangles all the files in the directory. It uses the - =config_tangle.el= file, which contains information required to - compute the current file names using for example ~(eval c)~ to get - the name of the produced C file. - - The file is not tangled if the last modification date of the org - file is less recent than one of the tangled files. - - #+begin_src bash -function tangle() -{ - local org_file=$1 - local c_file=${org_file%.org}.c - local f_file=${org_file%.org}.f90 - - if [[ ${org_file} -ot ${c_file} ]] ; then - return - elif [[ ${org_file} -ot ${f_file} ]] ; then - return - fi - emacs --batch ${org_file} --load=${top_srcdir}/tools/config_tangle.el -f org-babel-tangle -} - -for i in $@ -do - echo "--- ${i} ----" - tangle ${i} -done - #+end_src - -* Script to build the final qmckl.h file - :PROPERTIES: - :header-args:bash: :tangle build_qmckl_h.sh :noweb yes :shebang #!/bin/bash :comments org - :END: - - #+begin_src bash :noweb yes -# <> - - #+end_src - - #+NAME: qmckl-header - #+begin_src text :noweb yes ------------------------------------------- - QMCkl - Quantum Monte Carlo kernel library - ------------------------------------------ - - Documentation : <> - Issues : <> - - <> - - - #+end_src - - All the produced header files are concatenated in the =qmckl.h= - file, located in the include directory. The =*_private.h= files - are excluded. - - Put =.h= files in the correct order: - - #+begin_src bash -HEADERS="" -for i in $(cat table_of_contents) -do - HEADERS+="${i%.org}_type.h " -done - -for i in $(cat table_of_contents) -do - HEADERS+="${i%.org}_func.h " -done - #+end_src - - Generate C header file - - #+begin_src bash -OUTPUT="${top_srcdir}/include/qmckl.h" - -cat << EOF > ${OUTPUT} -/* - ,* <> - ,*/ - -#ifndef __QMCKL_H__ -#define __QMCKL_H__ - -#include -#include -#include -EOF - -for i in ${HEADERS} -do - if [[ -f $i ]] ; then - cat $i >> ${OUTPUT} - fi -done - -cat << EOF >> ${OUTPUT} -#endif -EOF - #+end_src - - Generate Fortran interface file from all =qmckl_*_fh.f90= files - - #+begin_src bash -HEADERS_TYPE="qmckl_*_fh_type.f90" -HEADERS="qmckl_*_fh_func.f90" - -OUTPUT="${top_srcdir}/share/qmckl/fortran/qmckl_f.f90" -cat << EOF > ${OUTPUT} -! -! <> -! -module qmckl - use, intrinsic :: iso_c_binding -EOF - -for i in ${HEADERS_TYPE} -do - cat $i >> ${OUTPUT} -done - -for i in ${HEADERS} -do - cat $i >> ${OUTPUT} -done - -cat << EOF >> ${OUTPUT} -end module qmckl -EOF - #+end_src - -* Script to build the documentation - :PROPERTIES: - :header-args:bash: :tangle build_doc.sh :noweb yes :shebang #!/bin/bash :comments org - :END: - - First define readonly global variables. - - #+begin_src bash :noweb yes -readonly DOCS=${top_srcdir}/share/doc/qmckl/ -readonly SRC=${top_srcdir}/src/ -readonly HTMLIZE=${DOCS}/html/htmlize.el -readonly CONFIG_DOC=${top_srcdir}/tools/config_doc.el -readonly CONFIG_TANGLE=${top_srcdir}/tools/config_tangle.el - #+end_src - - Check that all the defined global variables correspond to files. - - #+begin_src bash :noweb yes -function check_preconditions() -{ - if [[ -z ${top_srcdir} ]] - then - print "top_srcdir is not defined" - exit 1 - fi - - for dir in ${DOCS}/html ${DOCS}/text ${SRC} - do - if [[ ! -d ${dir} ]] - then - print "${dir} not found" - exit 2 - fi - done - - for file in ${CONFIG_DOC} ${CONFIG_TANGLE} - do - if [[ ! -f ${file} ]] - then - print "${file} not found" - exit 3 - fi - done -} - #+end_src - - ~install_htmlize~ installs the htmlize Emacs plugin if the - =htmlize.el= file is not present. - - #+begin_src bash :noweb yes -function install_htmlize() -{ - local url="https://github.com/hniksic/emacs-htmlize" - local repo="emacs-htmlize" - - [[ -f ${HTMLIZE} ]] || ( - cd ${DOCS}/html - git clone ${url} \ - && cp ${repo}/htmlize.el ${HTMLIZE} \ - && rm -rf ${repo} - cd - - ) - - # Assert htmlize is installed - [[ -f ${HTMLIZE} ]] \ - || exit 1 -} - #+end_src - - Extract documentation from an org-mode file. - - #+begin_src bash :noweb yes -function extract_doc() -{ - local org=$1 - local local_html=${SRC}/${org%.org}.html - local local_text=${SRC}/${org%.org}.txt - local html=${DOCS}/html/${org%.org}.html - local text=${DOCS}/text/${org%.org}.txt - - if [[ -f ${html} && ${org} -ot ${html} ]] - then - return - fi - emacs --batch \ - --load ${HTMLIZE} \ - --load ${CONFIG_DOC} \ - ${org} \ - --load ${CONFIG_TANGLE} \ - -f org-html-export-to-html \ - -f org-ascii-export-to-ascii - mv ${local_html} ${DOCS}/html - mv ${local_text} ${DOCS}/text - -} - #+end_src - - The main function of the script. - - #+begin_src bash :noweb yes -function main() { - - check_preconditions || exit 1 - - # Install htmlize if needed - install_htmlize || exit 2 - - # Create documentation - cd ${SRC} \ - || exit 3 - - for i in *.org - do - echo - echo "======= ${i} =======" - extract_doc ${i} - done - - if [[ $? -eq 0 ]] - then - cd ${DOCS}/html - rm -f index.html - ln README.html index.html - exit 0 - else - exit 3 - fi -} -main - #+end_src - - diff --git a/tools/build_makefile.py b/tools/build_makefile.py index 653497e..f85bcdb 100755 --- a/tools/build_makefile.py +++ b/tools/build_makefile.py @@ -48,7 +48,7 @@ def main(): c_test_o = "tests/test_"+i+".$(OBJEXT)" f_test_o = "tests/test_"+i+"_f.$(OBJEXT)" c_test = "tests/test_"+i+".c" - f_test = "tests/test_"+i+"_f.f90" + f_test = "tests/test_"+i+"_f.F90" html = "share/doc/qmckl/html/"+i+".html" text = "share/doc/qmckl/text/"+i+".txt" @@ -60,10 +60,10 @@ def main(): h_type=i+"_type.h" h_private_func=i+"_private_func.h" h_private_type=i+"_private_type.h" - f90=i+"_f.f90" + F90=i+"_f.F90" fo=i+"_f.$(OBJEXT)" - fh_func=i+"_fh_func.f90" - fh_type=i+"_fh_type.f90" + fh_func=i+"_fh_func.F90" + fh_type=i+"_fh_type.F90" ORG_FILES += [org] TANGLED_FILES += [tangled] @@ -132,17 +132,17 @@ def main(): DEPS[h_private_func] = [tangled] if "(eval f)" in grep: - F_FILES += [f90] + F_FILES += [F90] - if f90 in DEPS: - DEPS[f90] += [tangled, "$(src_qmckl_fo)"] + if F90 in DEPS: + DEPS[F90] += [tangled, "$(src_qmckl_fo)"] else: - DEPS[f90] = [tangled, "$(src_qmckl_fo)"] + DEPS[F90] = [tangled, "$(src_qmckl_fo)"] if fo in DEPS: - DEPS[fo] += [f90, "$(src_qmckl_fo)"] + DEPS[fo] += [F90, "$(src_qmckl_fo)"] else: - DEPS[fo] = [f90, "$(src_qmckl_fo)"] + DEPS[fo] = [F90, "$(src_qmckl_fo)"] if "(eval fh_func)" in grep: FH_FUNC_FILES += [fh_func] diff --git a/tools/build_qmckl_f.sh b/tools/build_qmckl_f.sh index 6ce896a..b58b023 100755 --- a/tools/build_qmckl_f.sh +++ b/tools/build_qmckl_f.sh @@ -1,9 +1,9 @@ #!/bin/sh -# Script to build the final src/qmckl_f.f90 file +# Script to build the final src/qmckl_f.F90 file set -e -# All the produced header files are concatenated in the =src/qmckl_f.f90= +# All the produced header files are concatenated in the =src/qmckl_f.F90= # file, located in the share/qmckl/fortran directory. @@ -30,8 +30,8 @@ fi # Generate Fortran interface file # ------------------------------- -HEADERS_TYPE="src/qmckl_*_fh_type.f90" -HEADERS="src/qmckl_*_fh_func.f90" +HEADERS_TYPE="src/qmckl_*_fh_type.F90" +HEADERS="src/qmckl_*_fh_func.F90" cat << EOF > ${src_qmckl_f} ! diff --git a/tools/config_tangle.el b/tools/config_tangle.el index 4c29fc9..8ce071d 100755 --- a/tools/config_tangle.el +++ b/tools/config_tangle.el @@ -39,15 +39,15 @@ (setq src (concat top_builddir "/src/")) (setq tests (concat top_builddir "/tests/")) (setq name (file-name-nondirectory (substring buffer-file-name 0 -4))) -(setq f (concat src name "_f.f90")) -(setq fh_func (concat src name "_fh_func.f90")) -(setq fh_type (concat src name "_fh_type.f90")) +(setq f (concat src name "_f.F90")) +(setq fh_func (concat src name "_fh_func.F90")) +(setq fh_type (concat src name "_fh_type.F90")) (setq c (concat src name ".c")) (setq h_func (concat src name "_func.h")) (setq h_type (concat src name "_type.h")) (setq h_private_type (concat src name "_private_type.h")) (setq h_private_func (concat src name "_private_func.h")) (setq c_test (concat tests "test_" name ".c")) -(setq f_test (concat tests "test_" name "_f.f90")) +(setq f_test (concat tests "test_" name "_f.F90")) (org-babel-lob-ingest (concat srcdir "/tools/lib.org")) diff --git a/tools/tangle.sh b/tools/tangle.sh index 2f37b95..1b1700b 100755 --- a/tools/tangle.sh +++ b/tools/tangle.sh @@ -22,7 +22,7 @@ function tangle() { local org_file=$1 local c_file=${org_file%.org}.c - local f_file=${org_file%.org}.f90 + local f_file=${org_file%.org}.F90 if [[ ${org_file} -ot ${c_file} ]] ; then return