Noncollinear Crystal Field reference frames

[MPerfetti]

Dear Easy spin community,
I am writing a program to calculate the magnetic properties of a single crystal composed by several (interacting) metal complexes.
The metals are lanthanides, therefore (in C1 symmetry) a full set of 27 CF parameters must be used. I have a set of parameters expressed in the xyz magnetic reference frame (the reference frame in which the ground state g tensor is diagonal). I know the Euler angles connecting this frame and the crystallographic reference frame (the frame in which the magnetic field is expressed and in which I want my output).
The two magnetic reference frames of the coupled lanthanides are noncollinear. I set up the spin system using Sys.B2Frame (B4Frame, B6Frame), however the program only runs if the frame contains 3 Euler angles, not 6 (i.e. the program rotates the CF frame but the frames of the two spins are always collinear).
If I do the same with a simple ZFS Hamiltonian (Sys.Dframe) everything works and Sys.Dframe must be given with 6 Euler angles.
I attach hereafter some examples (here I use only B20 for simplicity). I have also tried to set up B2 with 5 elements instead of 1 (i.e. [0 0 x 0 0] vs x) and it does not work.

VERSION 1_ TWO NONCOLLINEAR LN IONS UNCOUPLED. PROGRAM DOES NOT RUN
Sys.S=[15/2 15/2];
Sys.g=[4/3; 4/3];
Sys.J=0clight1e-4;
Sys.B2=[-1/3;-1/3]clight1e-4;
Sys.B2Frame=[0 0 0; 90 90 0]pi/180; %%with Sys.B2Frame=[0 0 0, 0 0 0] the program does not complain but the simulation is not correct
Bvec=[0 0 1];
Ham=ham(Sys,Bvec1000);

ERROR MESSAGE:
Operands to the logical AND (&&) and OR (||) operators must be convertible to logical scalar values.
Use the ANY or ALL functions to reduce operands to logical scalar values.
Error in validatespinsys
Error in ham
Error in Freeenergytorque (line 139)
Ham=ham(Sys,Bvec*1000);

VERSION 2_ TWO (COLLINEAR) LN IONS UNCOUPLED. PROGRAM RUNS BUT THE INFO ABOUT 2ND REF FRAME IS LOST, WRONG RESULT
Sys.S=[15/2 15/2];
Sys.g=[4/3; 4/3];
Sys.J=0clight1e-4;
Sys.B2=[-1/3;-1/3]clight1e-4;
Sys.B2Frame=[0 0 0]pi/180;
Bvec=[0 0 1];
Ham=ham(Sys,Bvec1000);

VERSION 3_TWO LN IONS COUPLED. PROGRAM WORKS AND THE SIMULATION IS CORRECT.
Sys.S=[15/2 15/2];
Sys.g=[4/3; 4/3];
Sys.J=0clight1e-4;
Sys.D=[-1;-1]clight1e-4;
Sys.DFrame=[0 0 0; 90 90 0]pi/180;
Bvec=[0 0 1];
Ham=ham(Sys,Bvec1000);

[Stefan Stoll]

This is a bug. Thanks for reporting.

[MPerfetti]

Dear Stefan,
Thanks for checking so fast.
I am trying to find a way to bypass the problem. A single spin seems to work (i.e. if I rotate the B2 ref frame the simulation is correct), so there must be a way to extract the rotated CF parameters. If I can access them, I can write a program that rotates the first ref frame, then rotates the second, and then construct a system with the rotated CF parameters.
Thank you,
Mauro

[Stefan Stoll]

It is now fixed in 6.0.1. Please have a look and post here to let us known whether it works as expected.

[MPerfetti]

Dear Stefan,
the new easy spin release (6.0.2) works fine. I have tested B2, B4 and B6 frame and they all seem to work. Thanks again for fixing the bug.
On this same issue, I wanted to ask if there is a way in easy spin to get:
1) the rotated CF parameters
2) the relation between the Euler angles of a certain spin and the others present in a certain space group (or equivalently, the Euler matrix of the other spins). At the moment the space group is only supported in the Exp structure, which I am not using in my program.

Thank you

[Stefan Stoll]

Unfortunately, this information is not directly accessible in EasySpin. You might want to look at the EasySpin source code on GitHub to find details about it. The relevant code is in ham_zf.m. Here is a code snipped:

Code: Select all

    % Apply transformation if non-zero tilt angles are given
    % (Sys.BFrame is processed from Sys.B?Frame by validatespinsys)
    tiltang = Sys.BFrame{k}(iSpin,:);
    if any(tiltang)
      % Calculate transformation matrix for ISTOs
      Dk = wignerd(k,tiltang(1),tiltang(2),tiltang(3));
      Ck = isto2stev(k);
      % Transform to Stevens operator basis
      DB = Ck*Dk.'/Ck;
      DB = removenumericalnoise(DB);
      % Transform Bk
      Bk_M = Bk(iSpin,:)*DB;  % eigenframe of Bk -> molecular frame
    else
      Bk_M = Bk(iSpin,:);
    end