!!------------------------------------------------------- !!---- Relaxed Selected Excitation (RelaxSE) !!------------------------------------------------------- !!---- This file is part of RelaxSE !!---- !!---- The RelaxSE project is distributed under LGPL. In agreement with the !!---- Intergovernmental Convention of the ILL, this software cannot be used !!---- in military applications. !!---- !!---- Copyright (C) 2016-2021 Institut Laue-Langevin (ILL), Grenoble, FRANCE !!---- Institut Neel - CNRS-UPR2940 (CNRS), Grenoble, FRANCE !!---- !!---- Authors: Elisa REBOLINI (ILL) rebolini@ill.fr !!---- Marie-Bernadette LEPETIT (CNRS) Marie-Bernadette.Lepetit@neel.cnrs.fr !!---- !!---- RelaxSE is free software; you can redistribute it and/or !!---- modify it under the terms of the GNU Lesser General Public !!---- License as published by the Free Software Foundation; either !!---- version 3.0 of the License, or (at your option) any later version. !!---- !!---- RelaxSE is distributed in the hope that it will be useful, !!---- but WITHOUT ANY WARRANTY; without even the implied warranty of !!---- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU !!---- Lesser General Public License for more details. !!---- !!---- You should have received a copy of the GNU Lesser General Public !!---- License along with this library; if not, see . !!---- !> @file Main for the properties of SASS program !> @author Marie-Bernadette Lepetit and Elisa Rebolini Program proprietes !!$ -------- Donness globales --------------------------------- use info use info_prop use files use typedet use spindetact use gener_ref0 use gener_ref1 use gener_monos use utils_char use utils_wrt #ifdef VAR_DEV use densite #endif !!$ -------- Donnes locales ----------------------------------- implicit none type(prop_infotype) :: prop_info type(prog_infotype) :: prog_info type(g_infotype) :: g_info type(o_infotype) :: o_info type(v_infotype) :: v_info type(det_infotype) :: det_info type(int_infotype) :: int_info type(sym_infotype) :: sym_info ! det type(deter), dimension(:), allocatable :: det, ref0 Integer (KIND=kd_int) :: ndet, nvec, nblock, nelact, nact Integer, dimension(:), allocatable :: shtblkdet, nblkdet, deter_index type(rlist) :: r type(spinrlist) :: rspin integer, parameter :: nb_rlist = 9 type(spindetact_list), pointer :: r0 type(deter_dblocklist) :: d ! vecteurs real(kd_dble), dimension(:,:), allocatable :: psi real(kd_dble), dimension(:), allocatable :: ener, psi_S ! prop Integer(KIND=kd_int) :: nprop character*5, dimension(:), allocatable :: whichprop ! proprietes real(kd_dble), dimension(:,:,:), allocatable :: rho_tot,rho_spin ! impressions integer, parameter :: pas = 10 character*1, parameter :: A1=" " ! autres integer :: idet, ivec, iprop, nb_thread integer :: nr, nq, iwr, i, j, nref0lu real(kd_dble) :: tstart, tend, wstart, wend, wstart1, wend1,t1, t2, wt1, wt2 CHARACTER(LEN=8) :: date ! returned values from DATE_AND_TIME() CHARACTER(LEN=10) :: time CHARACTER(LEN=5) :: zone INTEGER,DIMENSION(8) :: values !!$============================================================ !!$ -------- Code --------------------------------------------- !!$----- !!$----- Initialisations !!$----- call init_prop(g_info, prog_info, o_info, v_info, det_info, int_info, sym_info, prop_info) #ifdef VAR_OMP !$OMP PARALLEL SHARED(nb_thread) nb_thread=OMP_GET_MAX_THREADS() !$OMP END PARALLEL prog_info%nb_thread = nb_thread #endif ! Definition des fichiers call def_files_prop(prog_info, ndet) flush(f_output) call gettime(tstart, wstart) ! lecture des x_info read(f_info) prog_info read(f_info) g_info read(f_info) o_info read(f_info) v_info read(f_info) det_info call lire_intinfo(f_info, int_info) read(f_info) nblock allocate(shtblkdet(nblock), nblkdet(nblock), deter_index(nblock)) shtblkdet(:) = 0 nblkdet(:) = 0 deter_index = 99 read(f_info) shtblkdet(1:nblock), nblkdet(1:nblock), deter_index(1:nblock) read(f_info) sym_info call read_propinp(prop_info, nref0lu, f_input) prog_info%restart = .true. ndet = det_info%ndet nvec = v_info%nvec nelact = g_info%nelact nact = o_info%nact det_info%nref0 = nref0lu !!$----- !!$ --- Generation des det !!$----- if (prog_info%methodAct.eq.'cas') then det_info%nref0 = comb(o_info%nact,g_info%na) * comb(o_info%nact,g_info%nb) call deter_init(ref0, nref0lu) call gener_cas(ref0, prog_info, det_info, g_info, o_info) else call deter_init(ref0, nref0lu) call lect_ref0(ref0, prog_info, det_info, g_info, o_info, v_info) end if if (prog_info%id_cpu.eq.0) flush(f_output) call detact_all_init(r, nb_rlist) if (prog_info%id_cpu.eq.0) then write(f_output,*) write(f_output,*)">>> Generation of all determinants" flush(f_output) endif call gettime(t1,wstart1) ! Calcul des conf de spin des ref0 call generspin_ref0(ref0, vecref0, r, o_info, det_info, v_info, prog_info) ! Calcul des conf de spin des ref1 call generspin_ref1(ref0, r, o_info, det_info, g_info%nelact, prog_info) call deter_free(ref0) ! Generation of all active determinants call compute_all_spindetact(rspin, r, o_info, v_info%sz) ! Generation of all determinants call deter_blocklist_init(d) call generspin_det(det_info, rspin, d, o_info, v_info, prog_info, .false.) call gettime(t2,wend1) if (prog_info%id_cpu.eq.0) then write(f_output,*) '>>>Determinants generated in',t2-t1,'second(s)' flush(f_output) endif !!$----- !!$ --- Lecture des det !!$----- write(f_output,*) write(f_output,*) '>>> Reading determinants in file bdet' allocate (det(ndet)) call read_bdet(det,ndet,f_bdet) !!$----- !!$ --- Lecture des vecteurs !!$----- write(f_output,*) write(f_output,'(" >>> Reading vectors from file restart")') allocate(psi(ndet,nvec), ener(nvec)) psi(:,:) = 0.d0 ener(:) = 0.d0 call lect_vect(psi, ener, ndet, nvec) !!$----- !!$ --- Projection sur le cas !!$----- call gettime(t1,wt1) write(f_output,*) write(f_output,*) " >>> Vectors projection on active space" call cass_proj(psi, ener, ndet, nvec, nelact,det, det_info, o_info, & nblock, shtblkdet,nblkdet,deter_index,f_output) call gettime(t2,wt2) write(f_output,*) write(f_output,'(X,A,A,F17.2,A,F17.2,A)') '>>> Active space projection in ',& 'CPU time ',t2-t1,'s, Wall time ',wt2-wt1,' s' flush(f_output) !!$----- !!$ --- Autres propriétés !!$----- nprop = prop_info%nprop if (nprop.eq.0) goto 9999 allocate(whichprop(nprop)) whichprop(1:nprop) = prop_info%whichprop(1:nprop) do iprop = 1, nprop call lowercase(whichprop(iprop)) select case (whichprop(iprop)) case("s-s2") write(f_output,*) write(f_output,*) " >>> Total S " allocate(psi_S(nvec)) psi_S(:) = 0.d0 call s2(psi,psi_S,ndet,nvec,nelact, det, det_info, d, o_info, f_output) nq=nvec/10 nr=nvec - nq*pas do i =1,nq write(f_output,'(" S ")',advance='no') do iwr = 1, max(2*nact +1-6,0) write(f_output,'(a1)',advance='no') A1 end do write(f_output,9012) psi_S((i-1)*pas+1:i*pas) end do write(f_output,'(" S ")',advance='no') do iwr = 1, max(2*nact +1-6,0) write(f_output,'(a1)',advance='no') A1 end do write(f_output,9012) psi_S(nq*pas+1:nvec) write(f_output,*) deallocate(psi_S) 9012 format(10(F18.8,2x)) case("pref1") ! Projection sur les ref1 write(f_output,*) write(f_output,*) " >>> Projection on Ref1" write(f_output,*) " Not yet implemented" case("lcoef") ! Projection sur les ref1 write(f_output,*) write(f_output,*) " >>> Projection on determinants with coefficients largest than 0.05" write(f_output,*) " Not yet implemented" case("dens") ! Calcul de la matrice densite write(f_output,*) write(f_output,*) " >>> One particule density matrix" write(f_output,*) " Not yet implemented" #ifdef VAR_DEV call dens(ndet,nvec, psi,det, d,rspin, o_info,g_info, rho_tot,rho_spin) #endif end select end do 9001 format(a7," :",a80) 9993 format (5x,10(E15.6,1x)) 9994 format ("Ener ",10(F15.6,1x)) 9995 format (a5,10(F15.6,1x)) 9999 write(f_output,*) write(f_output,*) ' <<< End of prop code >>> ' write(f_output,*) deallocate(shtblkdet, nblkdet, deter_index) deallocate(det) deallocate(psi, ener) call gettime(tend,wend) call date_and_time(date, time, zone, values) write(f_output,'(X,A,F17.2,A,F17.2,A)') 'Calculation finished in CPUtime',& tend-tstart,'s Walltime:', & wend-wstart,'s' write(f_output,*) 'Calculation finished on ', date(7:8),'-',date(5:6),& '-',date(1:4), ' at ', time(1:2),':',time(3:4) end Program proprietes