- fixed a bug in option "-dislocation loop" where the loop was not co…

…rrectly constructed if it was not normal to the Z axis.

CHANGELOG

@@ -15,6 +15,7 @@ Master
 - fixed a bug where structure "C14" (with uppercase "C") was not recognized.
 - fixed a bug in mode "--polycrystal" that prevented constructing 2-D polycrystals normal to X or Y. Also fixed a bug where the seed was not sufficiently duplicated if the seed was not orthogonal; thanks to Zheyuan Xing, Southwest Jiaotong University, China, for pointing out this bug.
 - fixed a bug in option "-add-atom" that caused a segfault when the atom list was empty (not allocated).
+- fixed a bug in option "-dislocation loop" where the loop was not correctly constructed if it was not normal to the Z axis. Thanks to Marion Borde and Jonathan Amodeo, MATEIS, Univ. Lyon 1, France, for pointing out this bug.
 - fixed a bug in option "-select STL" when selecting atoms according to a complex 3-D shape from an STL file, where entire planes of atoms were not properly selected if the vertices were too close or if atom positions exactly matched the edges of several vertices.
 
 

src/modes/mode_polycrystal.f90

@@ -1076,25 +1076,49 @@ SUBROUTINE POLYCRYS(ucfile,vfile,options_array,prefix,outfileformats)
   templatebox(2) = templatebox(1)
   templatebox(3) = templatebox(1)
   templatebox(twodim) = VECLENGTH(H(twodim,:))
-ELSEIF( sameplane .OR. Nnodes<=4 ) THEN
-  !All nodes are in the same plane, or number of nodes is small
+ELSEIF( sameplane ) THEN
+  !All nodes are in the same plane
   !Increase template size to make sure all grains are covered
   templatebox(:) = 2.d0*templatebox(:)
+ELSEIF( Nnodes<=4 ) THEN
+  !Number of nodes is small
+  !Increase template size to make sure all grains are covered
+  templatebox(:) = 1.8d0*templatebox(:)
 ELSEIF( Nnodes>=10 ) THEN
   !Many nodes: decrease template size to improve performance
   templatebox(:) = templatebox(:) * MAX( 0.6d0 , 1.d0 - (Nnodes/200.d0) )
 ENDIF
 !Add a few angströms for good measure
 templatebox(:) = templatebox(:) + (/6.d0,6.d0,6.d0/)
+!Compute how many particles the template will contain = (seed density)*(template volume)
+P1 = CEILING( 1.1d0 * seed_density * PRODUCT(templatebox(:)) )
+WRITE(msg,'(a25,f18.0)') "Expected NP for template:", P1
+CALL ATOMSK_MSG(999,(/TRIM(msg)/),(/0.d0/))
+!If expected number of particles is too large, reduce template size
+IF( P1 > 2.147d9 ) THEN
+  templatebox(:) = 0.8d0*templatebox(:)
+  IF( twodim>0 ) THEN
+    templatebox(twodim) = VECLENGTH(H(twodim,:))
+  ENDIF
+ENDIF
+!Re-compute number of particles
+P1 = CEILING( 1.1d0 * seed_density * PRODUCT(templatebox(:)) )
+!If expected number of particles still is too large, abort completely
+IF( P1 > 2.147d9 ) THEN
+  CALL ATOMSK_MSG(821,(/""/),(/P1/))
+  nerr = nerr+1
+  GOTO 1000
+ENDIF
+!
 !By default the template grain is a bit larger than max.cell size * sqrt(3)
 !This template grain will be cut later to construct each grain
 expandmatrix(:) = 1
 WRITE(msg,*) "Determining expansion factors:"
 CALL ATOMSK_MSG(999,(/TRIM(msg)/),(/0.d0/))
 DO i=1,3
-  IF( VECLENGTH(Huc(i,:)) < 0.8d0*VECLENGTH(H(1,:)) .OR.  &
-    & VECLENGTH(Huc(i,:)) < 0.8d0*VECLENGTH(H(2,:)) .OR.  &
-    & VECLENGTH(Huc(i,:)) < 0.8d0*VECLENGTH(H(3,:))       ) THEN
+  IF( VECLENGTH(Huc(i,:)) < 0.8d0*MAXVAL(templatebox(:)) .OR.  &
+    & VECLENGTH(Huc(i,:)) < 0.8d0*MAXVAL(templatebox(:)) .OR.  &
+    & VECLENGTH(Huc(i,:)) < 0.8d0*MAXVAL(templatebox(:))       ) THEN
     !
     !Compute sum of seed vectors components along direction i
     P2 = DBLE( FLOOR( DABS( SUM(Huc(:,i)) )))
@@ -1112,44 +1136,14 @@ SUBROUTINE POLYCRYS(ucfile,vfile,options_array,prefix,outfileformats)
     expandmatrix(i) = CEILING(P1)
   ENDIF
 ENDDO
-WRITE(msg,*) "Initial expansion factors:", expandmatrix(:)
-CALL ATOMSK_MSG(999,(/TRIM(msg)/),(/0.d0/))
-!
 !If the system is 2-D, do not expand along the shortest axis
 IF( twodim>0 ) THEN
   expandmatrix(twodim) = 0
 ENDIF
-!Compute how many particles the template will contain
-!P1 = DBLE(MAX(1,expandmatrix(1))) * DBLE(MAX(1,expandmatrix(2))) * DBLE(MAX(1,expandmatrix(3))) * DBLE(SIZE(Puc,1))
-P1 = CEILING( 1.1d0 * seed_density * PRODUCT(templatebox(:)) )
-WRITE(msg,'(a25,f18.0)') "Expected NP for template:", P1
-CALL ATOMSK_MSG(999,(/TRIM(msg)/),(/0.d0/))
-!If expected number of particles is very large, reduce some values in expandmatrix(:)
-IF( P1 > 2.147d9 ) THEN
-  WRITE(msg,*) "NP is too large, reducing expansion factors...:"
-  CALL ATOMSK_MSG(999,(/TRIM(msg)/),(/0.d0/))
-  IF( Nnodes>10 ) THEN
-    !Many nodes in the box => reduce the size of template grain
-    expandmatrix(:) = NINT(0.7d0 * expandmatrix(:))
-  ELSE
-    !There are not many grains => do not reduce too much
-    expandmatrix(:) = NINT( 0.9d0 * expandmatrix(:) )
-  ENDIF
-ENDIF
-WRITE(msg,'(a47,3i6)') "Final corrected expansion factors for template:", expandmatrix(:)
+WRITE(msg,*) "Initial expansion factors:", expandmatrix(:)
 CALL ATOMSK_MSG(999,(/TRIM(msg)/),(/0.d0/))
-!Re-calculate expected number of atoms
-P1 = CEILING( 1.1d0 * seed_density * PRODUCT(templatebox(:)) )
-!P1 = DBLE(MAX(1,expandmatrix(1))) * DBLE(MAX(1,expandmatrix(2))) * DBLE(MAX(1,expandmatrix(3))) * DBLE(SIZE(Puc,1))
-!If expected NP is too large, abort completely
-IF( P1 > 2.147d9 ) THEN
-  CALL ATOMSK_MSG(821,(/""/),(/P1/))
-  nerr = nerr+1
-  GOTO 1000
-ENDIF
 !
 ! Allocate array newP for template grain (full duplicated crystal, not oriented or truncated)
-!m = PRODUCT(expandmatrix(:)+1)*SIZE(Puc,1)
 m = CEILING(P1)
 WRITE(msg,*) "ALLOCATE  newP, SIZE = ", m
 CALL ATOMSK_MSG(999,(/TRIM(msg)/),(/0.d0/))

src/options/disloc_loop.f90

@@ -254,4 +254,4 @@ END FUNCTION DisloSeg_displacement_iso
 !
 !
 !
-END MODULE dislocation_loop
+END MODULE dislocation_loop

src/options/opt_disloc.f90

@@ -1107,24 +1107,25 @@ SUBROUTINE DISLOC_XYZ(H,P,S,disloctype,dislocline,dislocplane,b,nu,pos,SELECT,OR
   !The direction normal to the loop is dislocline
   !pos(1:3) = position of the center of the loop; pos(4) = radius of the loop
   !
+  SELECT CASE(dislocline)
+  CASE('x','X')
+    !Loop in (y,z) plane
+    a1 = 2
+    a2 = 3
+    a3 = 1
+  CASE('y','Y')
+    a1 = 3
+    a2 = 1
+    a3 = 2
+  CASE DEFAULT
+    a1 = 1
+    a2 = 2
+    a3 = 3
+  END SELECT
+  !
   !Discretize loop into segments
   IF( pos(4) < 0.d0 ) THEN
     !"Negative" radius => user wants a square of side pos(4)
-    SELECT CASE(dislocline)
-    CASE('x','X')
-      !Loop in (y,z) plane
-      a1 = 2
-      a2 = 3
-      a3 = 1
-    CASE('y','Y')
-      a1 = 3
-      a2 = 1
-      a3 = 2
-    CASE DEFAULT
-      a1 = 1
-      a2 = 2
-      a3 = 3
-    END SELECT
     ALLOCATE(xLoop(4,3))
     xLoop(:,:) = 0.d0
     xLoop(1,a1) = pos(1) - pos(4)
@@ -1167,7 +1168,10 @@ SUBROUTINE DISLOC_XYZ(H,P,S,disloctype,dislocline,dislocplane,b,nu,pos,SELECT,OR
   ENDIF
   !For each atom, compute its displacement due to the loop
   DO i=1,SIZE(P,1)
-    P(i,1:3) = P(i,1:3) + LOOP_DISPLACEMENT(P(i,1:3), b, nu, pos(1:3), xLoop)
+    disp(:) = LOOP_DISPLACEMENT( P(i,1:3) , b, nu, pos(1:3) , xLoop)
+    P(i,a1) = P(i,a1) + disp(1)
+    P(i,a2) = P(i,a2) + disp(2)
+    P(i,a3) = P(i,a3) + disp(3)
   ENDDO
   !
   !