diff --git a/source/KKRnano/source/ProcessKKRresults_mod.F90 b/source/KKRnano/source/ProcessKKRresults_mod.F90
index 9e8f90f2bb6127376b3631831593e5d8d4439c2d..e003281c0999576bb9a7fb78f8f451e0edcfdfb8 100644
--- a/source/KKRnano/source/ProcessKKRresults_mod.F90
+++ b/source/KKRnano/source/ProcessKKRresults_mod.F90
@@ -1131,10 +1131,10 @@ module ProcessKKRresults_mod
reclen = max(reclen, 4*8 + 1*1 + 3*8)
end if
if (noncobfield) then
- ! The constraint bfields, torques, and moments before fixing directions
- ! (each double(3))
+ ! The constraint bfields (from this and from last iteration), torques,
+ ! and moments before fixing directions (each double(3))
recnum = recnum + 1
- reclen = max(reclen, (3+3+3)*8)
+ reclen = max(reclen, (3+3+3+3)*8)
end if
end subroutine
@@ -1206,8 +1206,11 @@ module ProcessKKRresults_mod
irec = irec + 1
endif
if (params%noncobfield) then
- write(unit=r1fu, rec=irec) calc%bfields(ila)%bfield_constr, &
- calc%bfields(ila)%mag_torque, calc%bfields(ila)%mag_mom
+ write(unit=r1fu, rec=irec) &
+ calc%bfields(ila)%bfield_constr, &
+ calc%bfields(ila)%last_bfield_constr, &
+ calc%bfields(ila)%mag_torque, &
+ calc%bfields(ila)%mag_mom
irec = irec + 1
end if
end do
@@ -1658,13 +1661,15 @@ module ProcessKKRresults_mod
double precision delta_angle !NOCO
integer :: max_delta_atom !NOCO
double precision, dimension(3) :: dir ! NOCO, direction of local frame of reference
- double precision, dimension(3) :: constr_field
+ double precision, dimension(3) :: constr_field, last_constr_field
double precision, dimension(3) :: torque ! torque of the current atom
double precision, dimension(3) :: newmoment ! magnetic moment before angle constraint
double precision :: max_torque_sqrd ! squared magnitude of largest torque so far
integer :: max_torque_atom ! atom with largest torque so far
double precision :: constr_angle_change, max_constr_angle_change
integer :: max_constr_angle_change_atom
+ double precision :: constr_field_change, max_constr_field_change
+ integer :: max_constr_field_change_atom
integer :: reclen, recnum, irec, ie, lrecres2
integer :: lcoremax, i1, ispin, lpot
character(len=*), parameter :: &
@@ -1680,6 +1685,8 @@ module ProcessKKRresults_mod
F82="(' Largest torque magnitude [Ry] = ',E12.6)", &
F83="(' Largest change of angle for constraint moment at atom = ', i5)", &
F84="(' Largest change of angle for constraint moment (deg) = ', f10.6)", &
+ F85="(' Largest change in constraint field at atom = ', i5)", &
+ F80="(' Largest change in constraint field (Ry) = ', e10.4)", &
F94="(4X,'nuclear charge ',F10.6,9X,'core charge = ',F10.6)"
integer :: npotd
@@ -1694,6 +1701,8 @@ module ProcessKKRresults_mod
max_torque_atom = 0
max_constr_angle_change = 0.
max_constr_angle_change_atom = -1 ! Initially a flag that means "feature not used (yet)"
+ max_constr_field_change = 0.
+ max_constr_field_change_atom = 0
! Calculate size of data written per atom
call calculateResults1FileShapes(reclen, recnum, compute_total_energy, &
@@ -1776,7 +1785,7 @@ module ProcessKKRresults_mod
irec = irec + 1
if (noncobfield) then
- read(71, rec=irec) constr_field, torque, newmoment
+ read(71, rec=irec) constr_field, last_constr_field, torque, newmoment
irec = irec + 1
end if
@@ -1831,6 +1840,17 @@ module ProcessKKRresults_mod
max_constr_angle_change = constr_angle_change
max_constr_angle_change_atom = i1
end if
+
+ ! Look for largest change in constraining field
+ ! Reuse variable for difference
+ last_constr_field = last_constr_field - constr_field
+ ! Calculate the absolute value of the change (squared)
+ constr_field_change = dot_product(last_constr_field, last_constr_field)
+ if (constr_field_change > max_constr_field_change) then
+ ! Update
+ max_constr_field_change = constr_field_change
+ max_constr_field_change_atom = i1
+ end if
end if
end if
@@ -1882,6 +1902,8 @@ module ProcessKKRresults_mod
if (max_constr_angle_change_atom >= 0) then
write(6, fmt=F83) max_constr_angle_change_atom
write(6, fmt=F84) 180.0/PI*max_constr_angle_change
+ write(6, fmt=F85) max_constr_field_change_atom
+ write(6, fmt=F80) sqrt(max_constr_field_change)
end if
end if
write(6, '(79(1h+))')
diff --git a/source/KKRnano/source/bfield/bfield.f90 b/source/KKRnano/source/bfield/bfield.f90
index 6b804b5a867b7612d3cf13c846cbeb24eeaf6c2a..c157fdfc2e661a7270bd74c75b8302e25bd7fffd 100644
--- a/source/KKRnano/source/bfield/bfield.f90
+++ b/source/KKRnano/source/bfield/bfield.f90
@@ -32,6 +32,7 @@ module mod_bfield
! Noncollinear magnetic fields
double precision, dimension(3) :: bfield_ext !! external magnetic field in cartesian coordinates
double precision, dimension(3) :: bfield_constr !! constraining field in cartesian coordinates
+ double precision, dimension(3) :: last_bfield_constr !! constraining field of last iteration
! Torque and other information used to update the constraints fields
double precision, dimension(3) :: mag_torque !! Magnetic torque
diff --git a/source/KKRnano/source/bfield/torque.f90 b/source/KKRnano/source/bfield/torque.f90
index a061c7ca94decb85388fe447f426ad2f4afd0872..51d3638f976ff35c64bb754e5fcb1586cfe08510 100644
--- a/source/KKRnano/source/bfield/torque.f90
+++ b/source/KKRnano/source/bfield/torque.f90
@@ -183,6 +183,9 @@ contains
mag_mom_dir = bfield%mag_mom(:) / mag_mom_len
mag_mom_sign = sign(1., dot_product(mag_mom_dir, dir))
+ ! Save previous iterations field
+ bfield%last_bfield_constr = bfield%bfield_constr
+
if (constr_mode == 3) then
bfield%bfield_constr(:) = bfield%bfield_constr(:) - &
( mag_mom_sign * (bfield%mag_torque(:) / mag_mom_len) * &