diff --git a/Code/Changes_private.txt b/Code/Changes_private.txt index 65a310cae4362040955fce0266e51e9565b7f8ef..b52840a033f8cf3ae9e17f4b106d1d672c637adc 100644 --- a/Code/Changes_private.txt +++ b/Code/Changes_private.txt @@ -8,3 +8,7 @@ for the private branch 'master' only: collected in 'user_PGR.f90'. These need explicit calls and coding to be activated. +- absorbing boundary conditions implemented + (special obseravbles as PES and PAD yet to be done). + +- protocol of energy and particle-number difference to initial state. diff --git a/Code/dynamic.f90 b/Code/dynamic.f90 index 915a7f28e1db8c76ae9294dc26fa6a46f8e41b3c..2897b93640b91e6012134e34566b3191709e40a4 100644 --- a/Code/dynamic.f90 +++ b/Code/dynamic.f90 @@ -24,6 +24,8 @@ MODULE DYNAMIC LOGICAL :: texternal=.FALSE. ! must be true if an external field is present LOGICAL :: text_timedep ! true for time-dependent external field REAL(db),PARAMETER :: esf=0.0D0 + REAL(db) :: pnrold(2)=(/0D0,0D0/),ehfold=0D0,eintold=0D0,ekinold=0D0,& + ecollold(2)=0D0 CONTAINS !************************************************************************* SUBROUTINE getin_dynamic @@ -49,7 +51,7 @@ CONTAINS ! Step 1: Preparation phase IF(.NOT.trestart) THEN iter=0 - time=0.0D0 + time=0D0 ! save wave functions CALL write_wavefunctions IF(wflag) WRITE(*,*) 'Wrote wave function file after initialization' @@ -65,6 +67,9 @@ CONTAINS CALL start_protocol(energiesfile, & '# Time N(n) N(p) E(sum) E(integ) Ekin & & Ecoll(n) Ecoll(p)') + CALL start_protocol(diffenergiesfile, & + '# Time N_n(0)-N_n(t) N_p(0)-N_p(t) diff-E(sum) diff-E(integ) diff-Ekin & + & diff-Ecoll(n) diff-Ecoll(p)') CALL start_protocol(monopolesfile, & '# Time rms_n rms_p rms_tot rms_n-rms_p') CALL start_protocol(dipolesfile, & @@ -151,7 +156,7 @@ CONTAINS ps4,rho,tau,current,sdens,sodens) psi(:,:,:,:,nst)=ps4 ENDDO - !$OMP END PARALLEL DO + !$OMP END PARALLEL DOx ! sum up over nodes IF(tmpi) CALL collect_densities IF(nabsorb > 0) CALL absbc(nabsorb,iter,nt,time) @@ -271,6 +276,13 @@ CONTAINS (current(:,:,:,1,iq)**2+current(:,:,:,2,iq)**2& +current(:,:,:,3,iq)**2)/rho(:,:,:,iq) ) ENDDO + IF(time==0D0) THEN + pnrold=pnr + ehfold=ehf + eintold=ehfint + ekinold=tke + ecollold=ecoll + END IF OPEN(unit=scratch,file=spinfile, POSITION='APPEND') WRITE(scratch,'(1x,i5,9F10.4)') iter,orbital,spin,total_angmom CLOSE(unit=scratch) @@ -278,6 +290,11 @@ CONTAINS WRITE(scratch,'(F10.2,2F8.3,2F15.7,3(1PG13.5))') & time,pnr,ehf,ehfint,tke,ecoll CLOSE(unit=scratch) + OPEN(unit=scratch,file=diffenergiesfile,POSITION='APPEND') + WRITE(scratch,'(F10.2,9(1pg13.5))') & + time,pnr-pnrold,ehf-ehfold,ehfint-eintold,& + tke-ekinold,ecoll-ecollold + CLOSE(unit=scratch) WRITE(*,'(/A,I7,A,F12.4,A/2(A,F12.4),A)') & ' ***** Time step #',iter,' at time ',time,' fm/c *****', & ' Total energy: ',ehf,' MeV Total kinetic energy: ', & diff --git a/Code/main3d.f90 b/Code/main3d.f90 index 15e0e9db753affa9ca7601499f586687956e9c96..5a56f767cc0858caba6635250eddb1486d8be2ae 100644 --- a/Code/main3d.f90 +++ b/Code/main3d.f90 @@ -16,7 +16,8 @@ PROGRAM tdhf3d INTEGER :: imode,nofsave !*********************************************************************** NAMELIST /files/ wffile,converfile,monopolesfile,dipolesfile, & - momentafile,energiesfile,quadrupolesfile,spinfile,extfieldfile + momentafile,energiesfile,quadrupolesfile,spinfile,extfieldfile,& + diffenergiesfile NAMELIST /main/ tcoul,mprint,mplot,trestart, & writeselect,write_isospin,mrest,imode,tfft,nof,r0 !******************************************************************** diff --git a/Code/params.f90 b/Code/params.f90 index 7a993e02dfd60cfa5fea48aaf341f03cd6465480..6a97da4e16cedde2328118acb15071bf4e242d2b 100644 --- a/Code/params.f90 +++ b/Code/params.f90 @@ -19,6 +19,7 @@ MODULE Params CHARACTER(LEN=80) :: dipolesfile='dipoles.res' CHARACTER(LEN=80) :: momentafile='momenta.res' CHARACTER(LEN=80) :: energiesfile='energies.res' + CHARACTER(LEN=80) :: diffenergiesfile='diffenergies.res' CHARACTER(LEN=80) :: quadrupolesfile='quadrupoles.res' CHARACTER(LEN=80) :: spinfile='spin.res' CHARACTER(LEN=80) :: extfieldfile='extfield.res'