Commit 12b9c6ea authored by Daniel Wortmann's avatar Daniel Wortmann

Added functionality to calculate E&B-Fields

parent 0cf0bbf4
......@@ -104,6 +104,7 @@ CONTAINS
REAL :: scale_temp, dtild_temp
REAL :: forceAllAtoms(3,atoms%nat)
CLASS(t_forcetheo),ALLOCATABLE:: forcetheo
TYPE(t_field) :: field
input=input_in
atoms_new=atoms
......@@ -207,7 +208,7 @@ CONTAINS
ALLOCATE(xmlCoreOccs(1,1,1))
CALL initWannierDefaults(wann_temp)
CALL r_inpXML(atoms_temp,obsolete_temp,vacuum_temp,input_temp,stars_temp,sliceplot_temp,&
banddos_temp,dimension_temp,forcetheo,cell_temp,sym_temp,xcpot_temp,noco_temp,&
banddos_temp,dimension_temp,forcetheo,field,cell_temp,sym_temp,xcpot_temp,noco_temp,&
oneD_temp,hybrid_temp,kpts_temp,enpara_temp,coreSpecInput_temp,wann_temp,noel_temp,&
namex_temp,relcor_temp,a1_temp,a2_temp,a3_temp,dtild_temp,xmlElectronStates,&
xmlPrintCoreStates,xmlCoreOccs,atomTypeSpecies,speciesRepAtomType,l_kpts_temp)
......
......@@ -258,54 +258,18 @@
INTEGER, INTENT(IN) :: k1d, k2d, nvac
REAL, INTENT(IN) :: area
REAL :: zsigma, sig_b(2), tmp
INTEGER, PARAMETER :: iou = 33
INTEGER :: ios
LOGICAL :: exists, eV
CHARACTER(len=50) :: str
INQUIRE (file='apwefl',exist=exists)
IF (.NOT. exists) RETURN
OPEN (iou,file='apwefl',form='formatted',status='old')
READ (iou, '(a)') str
str = ADJUSTL (str)
! Old format: zsigma, optionally followed by
! sig_b (1) and (2)
IF (str(1:1) /= '!' .AND. str(1:1) /= '#'&
& .AND. str(1:1) /= '&' .AND. str /= '') THEN
READ (str, *) E%zsigma
READ (iou, *, iostat=ios) sig_b
IF (ios == 0) E%sig_b = sig_b
RETURN
END IF
REAL :: tmp
INTEGER :: i
! New format
ALLOCATE(E%sigEF(3*k1d, 3*k2d, nvac))
E%sigEF = 0.0
eV = .false.
REWIND(iou)
DO
READ(iou,'(a)', IOSTAT=ios) str
IF (ios /= 0) EXIT
str = ADJUSTL (str)
IF (str(1:1) == '#' .OR. str(1:1) == '!' .OR. str == '')&
& CYCLE
IF (lower_case (str(1:8)) == '&efield ') THEN
BACKSPACE (iou)
CALL read_namelist (iou, E, eV)
ELSE IF (str(1:1) == '&') THEN
CALL juDFT_error("ERROR reading 'apwefl': Unknown namelist"&
& ,calledby ="efield")
ELSE
CALL read_shape (E, str, nvac)
END IF
DO i=1,SIZE(e%shapes)
CALL read_shape (E, e%shapes(i), nvac)
END DO
CLOSE (iou)
IF (eV) THEN
IF (e%l_eV) THEN
E%sig_b(:) = E%sig_b/hartree_to_ev_const
E%sigEF(:,:,:) = E%sigEF/hartree_to_ev_const
END IF
......
......@@ -14,7 +14,7 @@ MODULE m_rinpXML
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
CONTAINS
SUBROUTINE r_inpXML(&
atoms,obsolete,vacuum,input,stars,sliceplot,banddos,DIMENSION,forcetheo,&
atoms,obsolete,vacuum,input,stars,sliceplot,banddos,DIMENSION,forcetheo,field,&
cell,sym,xcpot,noco,oneD,hybrid,kpts,enpara,coreSpecInput,wann,&
noel,namex,relcor,a1,a2,a3,dtild,xmlElectronStates,&
xmlPrintCoreStates,xmlCoreOccs,atomTypeSpecies,speciesRepAtomType,&
......@@ -58,6 +58,7 @@ CONTAINS
TYPE(t_noco),INTENT(INOUT) :: noco
TYPE(t_dimension),INTENT(OUT) :: dimension
TYPE(t_enpara) ,INTENT(OUT) :: enpara
TYPE(t_field), INTENT(INOUT) :: field
CLASS(t_forcetheo),ALLOCATABLE,INTENT(OUT):: forcetheo
TYPE(t_coreSpecInput),INTENT(OUT) :: coreSpecInput
TYPE(t_wann) ,INTENT(INOUT) :: wann
......@@ -141,7 +142,7 @@ CONTAINS
INTEGER :: altKPointSetIndex, altKPointSetIndices(2)
LOGICAL :: ldaSpecies
REAL :: socscaleSpecies
REAL :: socscaleSpecies,b_field_mtspecies
INTEGER, ALLOCATABLE :: lNumbers(:), nNumbers(:), speciesLLO(:)
INTEGER, ALLOCATABLE :: loOrderList(:)
......@@ -700,27 +701,31 @@ input%preconditioning_param = evaluateFirstOnly(xmlGetAttributeValue('/fleurInpu
! END IF
! Read in optional E-Field parameters
xPathA = '/fleurInput/calculationSetup/eField'
xPathA = '/fleurInput/calculationSetup/fields'
numberNodes = xmlGetNumberOfNodes(xPathA)
IF (numberNodes.EQ.1) THEN
!input%efield%zsigma = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@zsigma'))
!input%efield%sig_b(1) = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@sig_b_1'))
!input%efield%sig_b(2) = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@sig_b_2'))
!input%efield%plot_charge = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@plot_charge'))
!input%efield%plot_rho = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@plot_rho'))
!input%efield%autocomp = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@autocomp'))
!input%efield%dirichlet = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@dirichlet'))
!l_eV = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@eV'))
STOP 'Error: Reading input for E-Fields not yet implemented completely!'
! ALLOCATE(input%efield%sigEF(3*k1d, 3*k2d, nvac))
! input%efield%sigEF = 0.0
!IF (l_eV) THEN
! input%efield%sig_b(:) = input%efield%sig_b/hartree_to_ev_const
! input%efield%sigEF(:,:,:) = input%efield%sigEF/hartree_to_ev_const
!END IF
IF (xmlGetNumberOfNodes(TRIM(ADJUSTL(xPathA))//'/@b_field')>0) THEN
field%b_field=evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'//@b_field'))
field%l_b_field=.true.
ENDIF
field%efield%zsigma = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@zsigma'))
field%efield%sig_b(1) = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@sig_b_1'))
field%efield%sig_b(2) = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@sig_b_2'))
field%efield%plot_charge = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@plot_charge'))
field%efield%plot_rho = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@plot_rho'))
field%efield%autocomp = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@autocomp'))
field%efield%dirichlet = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@dirichlet'))
field%efield%l_eV = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@eV'))
numberNodes=xmlGetNumberOfNodes(TRIM(ADJUSTL(xPathA))//'/shape')
ALLOCATE(field%efield%shapes(numberNodes))
DO i=1,numberNodes
WRITE(xPathB,"(a,a,i0,a)") TRIM(ADJUSTL(xpathA)),'/shape[',i,']'
field%efield%shapes(i)=TRIM(ADJUSTL(xmlGetAttributeValue(TRIM(ADJUSTL(xPathB)))))
ENDDO
END IF
! Read in optional energy parameter limits
......@@ -1409,6 +1414,10 @@ input%preconditioning_param = evaluateFirstOnly(xmlGetAttributeValue('/fleurInpu
IF (numberNodes==1) THEN
ldaSpecies = evaluateFirstBoolOnly(TRIM(ADJUSTL(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@lda'))))
socscaleSpecies = evaluateFirstOnly(TRIM(ADJUSTL(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@socscale'))))
IF (xmlGetNumberOfNodes(TRIM(ADJUSTL(xPathA))//'/@b_field_mt')>0) THEN
b_field_mtSpecies=evaluateFirstOnly(TRIM(ADJUSTL(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@b_field_mt'))))
field%l_b_field=.TRUE.
ENDIF
ENDIF
! Explicitely provided core configurations
......@@ -1621,6 +1630,13 @@ input%preconditioning_param = evaluateFirstOnly(xmlGetAttributeValue('/fleurInpu
xcpot%lda_atom(iType)=ldaSpecies
END SELECT
noco%socscale(iType)=socscaleSpecies
IF (field%l_b_field) THEN
IF (.NOT.ALLOCATED(field%b_field_mt)) THEN
ALLOCATE(field%b_field_mt(atoms%ntype))
field%b_field_mt=0.0
ENDIF
field%b_field_mt(itype)=b_field_mtSpecies
ENDIF
END IF
END DO
DEALLOCATE(loOrderList)
......
......@@ -52,7 +52,7 @@
<xsd:element maxOccurs="1" minOccurs="0" name="geometryOptimization" type="GeometryOptimizerType"/>
<xsd:element maxOccurs="1" minOccurs="0" name="ldaU" type="LdaUGeneralType"/>
<xsd:element maxOccurs="1" minOccurs="0" name="spinSpiralQPointMesh" type="SpinSpiralQPointMeshType"/>
<xsd:element maxOccurs="1" minOccurs="0" name="eField" type="EFieldType"/>
<xsd:element maxOccurs="1" minOccurs="0" name="fields" type="FieldType"/>
<xsd:element maxOccurs="1" minOccurs="0" name="energyParameterLimits" type="EnergyParameterLimitsType"/>
</xsd:all>
</xsd:complexType>
......@@ -62,7 +62,11 @@
<xsd:attribute name="elup" type="xsd:string"/>
</xsd:complexType>
<xsd:complexType name="EFieldType">
<xsd:complexType name="FieldType">
<xsd:sequence>
<xsd:element maxOccurs="unbounded" minOccurs="0" name="shape" type="xsd:string" />
</xsd:sequence>
<xsd:attribute default="0.0" name="b_field" type="xsd:string" use="optional"/>
<xsd:attribute default="10.0" name="zsigma" type="xsd:string" use="optional"/>
<xsd:attribute default="0.0" name="sig_b_1" type="xsd:string" use="optional"/>
<xsd:attribute default="0.0" name="sig_b_2" type="xsd:string" use="optional"/>
......@@ -356,6 +360,8 @@
<xsd:complexType name="SpecialType">
<xsd:attribute name="lda" type="FleurBool" use="optional" default="F"/>
<xsd:attribute name="socscale" type="ZeroToOneNumberType" use="optional" default="1.0"/>
<xsd:attribute name="b_field_mt" type="xsd:string" use="optional" default="0.0"/>
<xsd:attribute name="vca_charge" type="xsd:string" use="optional" default="0.0"/>
</xsd:complexType>
......
......@@ -172,7 +172,7 @@
a2 = 0.0
a3 = 0.0
CALL r_inpXML(&
atoms,obsolete,vacuum,input,stars,sliceplot,banddos,DIMENSION,forcetheo,&
atoms,obsolete,vacuum,input,stars,sliceplot,banddos,DIMENSION,forcetheo,field,&
cell,sym,xcpot,noco,oneD,hybrid,kpts,enpara,coreSpecInput,wann,&
noel,namex,relcor,a1,a2,a3,dtild,xmlElectronStates,&
xmlPrintCoreStates,xmlCoreOccs,atomTypeSpecies,speciesRepAtomType,&
......
......@@ -24,6 +24,8 @@ MODULE m_types_field
LOGICAL :: autocomp = .TRUE. ! Auto-compensate film charge
LOGICAL :: dirichlet = .FALSE. ! Dirichlet vs. Neumann boundary cond.
LOGICAL :: l_dirichlet_coeff = .FALSE. ! For MPI, true if C1/C2 set
LOGICAL :: l_eV =.FALSE. !Input in eV
CHARACTER(len=50),ALLOCATABLE :: shapes(:)
END TYPE t_efield
TYPE t_field
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment