Rearranged commit ee1d1e88.

This looks a little bit nicer.

inpgen/crystal.f

@@ -60,8 +60,8 @@
       INTEGER, PARAMETER  :: neig12=12
 
 !===> Local Variables
-      LOGICAL lerr,err_setup,invsym
-      INTEGER i,j,n,m,na,nt,mdet,mtr,nop0,fh
+      LOGICAL lerr,err_setup,invsym,inversionOp
+      INTEGER i,j,k,n,m,na,nt,mdet,mtr,nop0,fh
       REAL    t,volume,eps7,eps12
 
       INTEGER optype(nop48)
@@ -170,6 +170,44 @@
      >                atomid,atompos,natin,nop48,neig12,
      <                ntype,nat,nops,mrot,tau,
      <                neq,ntyrep,zatom,natype,natrep,natmap,pos)
+         ! Check whether there is an inversion center that is not at the
+         ! origin and if one is found shift the crystal such that the
+         ! inversion is with respect to the origin. Then recalculate
+         ! symmetry operations.
+         inversionOp = -1
+         symOpLoop: DO k = 1, nops
+            DO i = 1, 3
+               DO j = 1, 3
+                  IF (i.EQ.j) THEN
+                     IF (mrot(i,j,k).NE.-1) CYCLE symOpLoop
+                  ELSE
+                     IF (mrot(i,j,k).NE.0) CYCLE symOpLoop
+                  END IF
+                  IF ((i.EQ.3).AND.(j.EQ.3)) THEN
+                     inversionOp = k
+                     EXIT symOpLoop
+                  END IF
+               END DO
+            END DO
+         END DO symOpLoop
+         IF (inversionOp.GT.0) THEN
+            IF(ANY(ABS(tau(:,inversionOp)).GT.eps7)) THEN
+               WRITE(*,*) 'Found inversion center at finite position.'
+               WRITE(*,*) 'Shifting crystal by:'
+               WRITE(*,'(3f15.10)') 0.5*tau(:,inversionOp)
+               WRITE(*,*) ''
+               DO k = 1, ABS(natin)
+                  atompos(:,k) = atompos(:,k) + 0.5*tau(:,inversionOp)
+               END DO
+               DEALLOCATE(neq,ntyrep,zatom,mrot,tau)
+               CALL spg_gen(
+     >                      dispfh,outfh,errfh,dispfn,
+     >                      .FALSE.,symor,as,bs,scale,
+     >                      atomid,atompos,natin,nop48,neig12,
+     <                      ntype,nat,nops,mrot,tau,
+     <                      neq,ntyrep,zatom,natype,natrep,natmap,pos)
+            END IF
+         END IF
       ENDIF
 
       WHERE ( abs( tau ) < eps7 ) tau = 0.00

inpgen/spg_gen.f

@@ -43,7 +43,7 @@
       INTEGER index_op(nop48),num_tr(nop48)
 
       INTEGER mtmp(3,3),mp(3,3)
-      INTEGER u,v,w, inversionOp
+      INTEGER u,v,w
       INTEGER i,j,k,n,mop,nc,new,ns,nt,ntypm,ind(1),iTr,maxTrVecs
       INTEGER ity(natin),npaint,ipaint(natin)
       INTEGER ios,istep0, binSize, maxBinSize
@@ -70,22 +70,6 @@
      >                  as,bs,scale,nop48,neig12,
      <                  mops,mmrot)
 
-!---> read in atomic positions and shift to (-1/2,1/2] in lattice
-!---> coords. also read in identification (atomic) number (atomid)
-!---> to distinguish different atom types (need not be atomic number)
-      DO n=1,natin
-         IF (cartesian) THEN  ! convert to lattice coords. if necessary
-            atompos(:,n) = matmul( bs, atompos(:,n) )
-         ENDIF
-      END DO
-
-  200 CONTINUE !The routine jumps back to this point with shifted atom positions
-               !an inversion center is found that is not the origin
-
-      DO n=1,natin
-         pos(:,n) = atompos(:,n) - anint( atompos(:,n) - eps7 )
-      END DO
-
       ALLOCATE ( mtable(mops,mops) )
       CALL close_pt(                        ! check closure and get
      >              mops,mmrot,             ! multiplication table
@@ -99,6 +83,16 @@
          ENDDO
       ENDDO
 
+!---> read in atomic positions and shift to (-1/2,1/2] in lattice
+!---> coords. also read in identification (atomic) number (atomid)
+!---> to distinguish different atom types (need not be atomic number)
+      DO n=1,natin
+         IF (cartesian) THEN  ! convert to lattice coords. if necessary
+            atompos(:,n) = matmul( bs, atompos(:,n) )
+         ENDIF
+         pos(:,n) = atompos(:,n) - anint( atompos(:,n) - eps7 )
+      ENDDO
+
 !---> store the positions (in lattice coord.s) given in the input file
 !     ALLOCATE ( inipos(3,natin) )
 !     DO n=1,natin
@@ -514,8 +508,6 @@
 !      if( nopd < nops )then
 !        nopd = nops
 !      endif
-      IF (ALLOCATED(mrot)) DEALLOCATE(mrot)
-      IF (ALLOCATED(tau)) DEALLOCATE(tau)
       ALLOCATE ( mrot(3,3,nops),tau(3,nops) )
 
       DO n=1,nops
@@ -572,42 +564,9 @@
          ENDDO
       ENDDO
 
-      ! Check whether there is an inversion center that is not the origin
-      ! and shift the crystal if there is one. Then calculate symmetries again.
-      inversionOp = -1
-      symOpLoop: DO k = 1, nops
-         DO i = 1, 3
-            DO j = 1, 3
-               IF (i.EQ.j) THEN
-                  IF (mrot(i,j,k).NE.-1) CYCLE symOpLoop
-               ELSE
-                  IF (mrot(i,j,k).NE.0) CYCLE symOpLoop
-               END IF
-               IF ((i.EQ.3).AND.(j.EQ.3)) THEN
-                  inversionOp = k
-                  EXIT symOpLoop
-               END IF
-            END DO
-         END DO
-      END DO symOpLoop
-      IF (inversionOp.GT.0) THEN
-         IF(ANY(ABS(tau(:,inversionOp)).GT.eps7)) THEN
-            WRITE(*,*) 'Found inversion center at finite position.'
-            WRITE(*,*) 'Shifting crystal by:'
-            WRITE(*,'(3f15.10)') 0.5*tau(:,inversionOp)
-            WRITE(*,*) ''
-            DO k = 1, ABS(natin)
-               atompos(:,k) = atompos(:,k) + 0.5*tau(:,inversionOp)
-            END DO
-            GO TO 200
-         END IF
-      END IF
-
-
 !      if( ntypd < ntype )then
 !        ntypd = ntype
 !      endif
-
       ALLOCATE( neq(ntype),ntyrep(ntype),zatom(ntype) )
 
       neq(1:ntype) = 0