rw_noco.f90 5.74 KB
Newer Older
1 2 3 4 5 6
      MODULE m_rwnoco
      USE m_juDFT
!---------------------------------------------------------------------
!     read or write the nocoinp-file
!---------------------------------------------------------------------
      CONTAINS
7
        SUBROUTINE rw_noco_read(atoms,noco,input)
8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
          USE m_types
          USE m_constants
          IMPLICIT NONE
          ! ..
          ! ..  Arguments ..
          TYPE(t_atoms),INTENT(IN)   :: atoms
          TYPE(t_noco),INTENT(INOUT) :: noco
          TYPE(t_input),INTENT(INOUT):: input     


          ! ..
          INTEGER :: itype, j, fileend
          REAL    :: qsc(3)
          CHARACTER(len=8) :: inpchar 



          DO itype = 1,atoms%ntype
26 27
             READ (24,'(22x,l1)') noco%l_relax(itype)
             
28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63
             inpchar(1:2)= 'XX'
             READ (24,fmt='(19x,a2)') inpchar(1:2)
             BACKSPACE (24)
             IF (inpchar(1:2)/='pi') THEN
                READ (24,fmt='(7x,f14.10,11x,f14.10)') &
                     &      noco%alph(itype),noco%b_con(1,itype)
             ELSE
                READ (24,fmt='(7x,f12.8,13x,f14.10)') &
                     &      noco%alph(itype),noco%b_con(1,itype)
                noco%alph(itype)= noco%alph(itype)*pi_const
             ENDIF
             inpchar(1:2)= 'XX'
             READ (24,fmt='(19x,a2)') inpchar(1:2)
             BACKSPACE (24)
             IF (inpchar(1:2)/='pi') THEN
                READ (24,fmt='(7x,f14.10,11x,f14.10)') &
                     &      noco%beta(itype),noco%b_con(2,itype)
             ELSE
                READ (24,fmt='(7x,f12.8,13x,f14.10)') &
                     &      noco%beta(itype),noco%b_con(2,itype)
                noco%beta(itype)= noco%beta(itype)*pi_const
             ENDIF

             READ (24,*)
                WRITE (6,8010)        itype,noco%l_relax(itype)
                WRITE (6,8020)        noco%alph(itype),noco%b_con(1,itype)
                WRITE (6,8025)         noco%beta(itype),noco%b_con(2,itype)
                WRITE (6,*)
          ENDDO
8010      FORMAT ('atom-type',i4,',l_relax=',l1,',l_magn=',l1,&
               &',M=',f8.5,',magtype=',i4)
8020      FORMAT ('alpha =',f14.10,',b_cons_x =',f14.10)
8025      FORMAT ('beta  =',f14.10,',b_cons_y =',f14.10)
8026      FORMAT ('Total number of magnetic atoms',i4,',magnetic types',i4)

          READ (24,*)
64
          READ (24,8036) noco%l_ss,noco%l_mperp,noco%l_constr
65 66 67 68 69 70 71
!!$          IF ( (inpchar(1:8)=='sso_opt=') .OR. (noco%l_ss .AND. noco%l_soc) ) THEN
!!$             BACKSPACE (24)
!!$             READ (24,fmt='(45x,2l1)') input%sso_opt(1),input%sso_opt(2)
!!$          ELSE
!!$             input%sso_opt(1)= .FALSE. 
!!$             input%sso_opt(2)= .FALSE.
!!$          ENDIF
72
          READ (24,8046) noco%mix_b
73 74

          WRITE (6,fmt='(5(A,l1),l1)') &
75 76
               & 'l_ss=',noco%l_ss,',l_mperp=',noco%l_mperp,',l_constr=',noco%l_constr
          WRITE (6,8040) noco%mix_b
77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111
8030      FORMAT ('l_ss=',l1,',l_mperp=',l1,',l_constr=',l1,',l_disp=',l1)
8035      FORMAT (5x,l1,9x,l1,10x,l1,8x,l1)
8036      FORMAT (5x,l1,9x,l1,10x,l1)
8040      FORMAT ('mix_b=',f6.3,',thetaJ=',f14.10,',nsh=',i4)
8045      FORMAT (6x,f6.3,8x,f14.10,5x,i4)
8046      FORMAT (6x,f6.3)

          IF (noco%l_ss) THEN
             READ (24,fmt='(5x,3(f14.10,1x))') &
                  &    noco%qss(1),noco%qss(2),noco%qss(3)
             inpchar(1:3)= 'XXX'
             READ (24,fmt='(a4)',iostat=fileend) inpchar(1:4)
             IF (fileend==0) THEN
                IF (inpchar(1:4)=='qsc=') THEN
                   BACKSPACE (24)
                   READ (24,fmt='(5x,3(f14.10,1x))') &
                        &      qsc(1),qsc(2),qsc(3)
                   DO j= 1,3
                      IF ( ABS(qsc(j)) < 1.e-6 ) THEN
                         WRITE (6,fmt='(A,i1,A,1x,f14.10)')&
                              &          'Error reading nocoinp: qsc(',j,') =',qsc(j)
                         CALL juDFT_error("Error reading nocoinp",calledby&
                              &              ="rw_noco")
                      ENDIF
                      noco%qss(j)= noco%qss(j)/qsc(j)
                   ENDDO
                ENDIF
             ENDIF
             WRITE (6,*) 'This is a Spin-Spiral (SS) calculation. The'
             WRITE (6,*) 'q-vector of the Spin-Spiral is:'
             WRITE (6,8060) noco%qss(1),noco%qss(2),noco%qss(3)
8060         FORMAT ('qss=(',f14.10,',',f14.10,',',f14.10,')')
          ENDIF
        END SUBROUTINE rw_noco_read

112
      SUBROUTINE rw_noco_write(atoms,noco,input)
113 114 115 116 117 118 119 120 121 122 123 124 125 126
      USE m_types
      IMPLICIT NONE
! ..
! ..  Arguments ..
      TYPE(t_atoms),INTENT(IN)   :: atoms
      TYPE(t_noco),INTENT(INOUT) :: noco
      TYPE(t_input),INTENT(IN)   :: input     

! ..
      INTEGER :: itype
      CHARACTER(len=8) :: inpchar 
 
      DO itype = 1,atoms%ntype
         IF ( noco%l_relax(itype) ) noco%b_con(:,itype) = 0.0
127
         WRITE (24,8010) itype,noco%l_relax(itype)
128 129 130 131 132 133 134
         WRITE (24,8020) noco%alph(itype),noco%b_con(1,itype)
         WRITE (24,8025)  noco%beta(itype),noco%b_con(2,itype)
         WRITE (24,*)
      ENDDO
 
      WRITE (24,*) '-- logical parameters --'
      WRITE (24,fmt='(5(A,l1),l1)') &
135 136
     & 'l_ss=',noco%l_ss,',l_mperp=',noco%l_mperp,',l_constr=',noco%l_constr
      WRITE (24,8040) noco%mix_b
137 138 139 140 141 142 143 144 145 146 147 148 149

      IF (noco%l_ss) THEN
       WRITE (24,8060) noco%qss(1),noco%qss(2),noco%qss(3)
       WRITE (6,8060) noco%qss(1),noco%qss(2),noco%qss(3)
      ENDIF
8010      FORMAT ('atom-type',i4,',l_relax=',l1,',l_magn=',l1,&
               &',M=',f8.5,',magtype=',i4)
8020      FORMAT ('alpha =',f14.10,',b_cons_x =',f14.10)
8025      FORMAT ('beta  =',f14.10,',b_cons_y =',f14.10)
8040      FORMAT ('mix_b=',f6.3,',thetaJ=',f14.10,',nsh=',i4)
8060         FORMAT ('qss=(',f14.10,',',f14.10,',',f14.10,')')
      END SUBROUTINE rw_noco_write
      END MODULE m_rwnoco