Cartoon of spin-orbit coupling 

Matlab code 

% Cartoon of the spin-orbit interaction


clear

close all

vid = VideoWriter('SpinOrbitCoupling.mp4','MPEG-4');

vid.Quality = 100;

vid.FrameRate = 30;

open(vid);

figure('units','pixels','position',[0 0 1920 1080],'ToolBar','none');

set(0,'defaultfigurecolor',[1 1 1]);

set(gca,'linewidth',7); 


% Plotting areas of two components of video 

pos1 = [0 -0.1 0.7 1.28]; % Cartoon electron orbiting nucleus 

pos2 = [0.7 0.25 0.25 0.5]; % Energy-level scheme 


lp = [7 5 8]; % Angle of lighting

myNickel = [0.44 0.46 0.44];

myBlue = [0.4 0.44 0.73]; 


sRotAlpRed(1,:) = 0:2:3*360-2; % Rotation 2-degree steps over three rotations 

alp1 = zeros(1,90); % No central magnet for first half rotation (90 frames)

alp2 = 0:0.05:1; % Central magnet fade in over 21 frames 

alp3 = ones(1,3*180-90-21); % Central magnet 100% visible

sRotAlpRed(2,:) = [alp1 alp2 alp3]; % Second row determines behaviour of transparency of central magnet 

alp1 = zeros(1,180); % No spin magnet for first rotation (180 frames)

alp2 = 0:0.05:1; % Spin magnet fade in over 21 frames 

alp3 = ones(1,3*180-180-21); % Spin magnet 100% visible

sRotAlpRed(3,:) = [alp1 alp2 alp3]; % Third row determines behaviour of transparency of spin magnet 

Red1 = zeros(1,360); % No magnets for first rotation then in second rotation the colour is flipped from red to black

Red2 = 0:0.05:1; % s-magnet colour changes from black to red over 21 frames

Red3 = ones(1,3*180-360-21); 

sRotAlpRed(4,:) = [Red1 Red2 Red3]; % Fourth row determines colour flip for spin magnet 


[a,b,c] = sphere(50); % Electron or nucleon


rA = 0:0.125:0.125; % For cone arrowhead

[X,Y,Z] = cylinder(rA,50); % Arrow head

arrowLength = 0.4; % Length of arrow head


r2 = 0.22; % Radius of electron 

fac = 1.4; % Field-line ellipse b/a ratio 

gr = (1 + sqrt(5))/2; % Golden ratio 

nuclCol = [0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1]; % Red or black nucleons 

magFalp = 0.5; % Transparency of magnetic field lines 


R=5; % outer radius of torus representinig electron orbit 


for elRot = sRotAlpRed(1,:) % Three rotations 

    index = (elRot+2)/2; 

    subplot('position',pos1); % Cartoon of electron orbiting nucleus 

    newplot

    nucIndex = 0;

    hold off

    % Create filled dodecahedron defined by spheres at the vertices for

    % nuclear core 

    r = 0.16;

    for i = (-1:2:1)*r*fac

        for j = (-1:2:1)*r*fac

            for k = (-1:2:1)*r*fac

                nucIndex = nucIndex+1;

                if (nuclCol(nucIndex)>=0.5)

                    borr = 'r';

                else

                    borr = 'k';

                end

                surf(r*a+i,r*b+j,r*c+k,...

                    'FaceColor',borr,'LineStyle','none','FaceAlpha',1-sRotAlpRed(2,index),...

                    'FaceLighting','flat','DiffuseStrength',1);

                if (nuclCol(nucIndex)>=0.5)

                    borr = 'r';

                else

                    borr = 'k';

                end

                surf(r*a+i/fac,r*b+j/fac,r*c+k/fac,...

                    'FaceColor','r','LineStyle','none','FaceAlpha',(1-sRotAlpRed(2,index)),...

                    'FaceLighting','flat','DiffuseStrength',1);

                hold on

            end

        end

    end

    for i = (-gr:2*gr:gr)*r*fac

        for j = (-1/gr:2/gr:1/gr)*r*fac

            nucIndex = nucIndex+1;

            if (nuclCol(nucIndex)>=0.5)

                borr = 'r';

            else

                borr = 'k';

            end

            surf(r*a,r*b+i,r*c+j,...

                'FaceColor',borr,'LineStyle','none','FaceAlpha',(1-sRotAlpRed(2,index)),...

                'FaceLighting','flat','DiffuseStrength',1);

            nucIndex = nucIndex+1;

            if (nuclCol(nucIndex)>=0.5)

                borr = 'r';

            else

                borr = 'k';

            end

            surf(r*a+j,r*b,r*c+i,...

                'FaceColor',borr,'LineStyle','none','FaceAlpha',(1-sRotAlpRed(2,index)),...

                'FaceLighting','flat','DiffuseStrength',1);

            nucIndex = nucIndex+1;

            if (nuclCol(nucIndex)>=0.5)

                borr = 'r';

            else

                borr = 'k';

            end

            surf(r*a+i,r*b+j,r*c,...

                'FaceColor',borr,'LineStyle','none','FaceAlpha',(1-sRotAlpRed(2,index)),...

                'FaceLighting','flat','DiffuseStrength',1);

            hold on

        end

    end

    

    % Make central magnet representing angular momentum l 

    magWidth = 0.7; 

    [V,F] = platonic_solid(2,magWidth); % Octahedron, size = 1

    V(:,3) = sqrt(3)*V(:,3)-magWidth;

    patch('Faces',F,'Vertices',V,'FaceColor','r','FaceAlpha',sRotAlpRed(2,index), ...

        'EdgeColor','none','FaceLighting','flat','DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    [V,F] = platonic_solid(2,magWidth); % Octahedron, size = 1

    V(:,3) = sqrt(3)*V(:,3)+magWidth;

    patch('Faces',F,'Vertices',V,'FaceColor','k','FaceAlpha',sRotAlpRed(2,index), ...

        'EdgeColor','none','FaceLighting','flat','DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    material shiny


    % Make electron magnet representing angular momentum s

    magWidth = 0.4; 

    [V,F] = platonic_solid(2,magWidth); % Octahedron, size = 1

    V(:,1) = V(:,1)-R;

    V(:,3) = sqrt(3)*V(:,3) + magWidth;

    colNow = [1-sRotAlpRed(4,index) 0 0]; 

    ps = patch('Faces',F,'Vertices',V,'FaceColor',colNow,'FaceAlpha',sRotAlpRed(3,index), ...

        'EdgeColor','none','FaceLighting','flat','DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    direction = [0 0 1];

    rotate(ps,direction,-elRot,[-R 0 0])

    rotate(ps,direction,elRot,[0 0 0])

    [V,F] = platonic_solid(2,magWidth); % Octahedron, size = 1

    V(:,1) = V(:,1)-R;

    V(:,3) = sqrt(3)*V(:,3)-magWidth;

    colNow = [sRotAlpRed(4,index) 0 0]; 

    ps = patch('Faces',F,'Vertices',V,'FaceColor',colNow,'FaceAlpha',sRotAlpRed(3,index), ...

        'EdgeColor','none','FaceLighting','flat','DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    direction = [0 0 1];

    rotate(ps,direction,-elRot,[-R 0 0])

    rotate(ps,direction,elRot,[0 0 0])

    material shiny


    

    % Create a very thin torus to describe the electron orbit

    r=0.032; % inner tube radius

    th=linspace(0,2*pi,90); % e.g. 90 partitions along perimeter of the tube

    phi=linspace(0,2*pi,180); % e.g. 180 partitions along azimuth of torus

    % we convert our vectors phi and th to [n x n] matrices with meshgrid command:

    [Phi,Th]=meshgrid(phi,th);

    % now we generate n x n matrices for x,y,z according to eqn of torus

    x=(R+r.*cos(Th)).*cos(Phi);

    y=(R+r.*cos(Th)).*sin(Phi);

    z=r.*sin(Th);

    surf(x,y,z,'FaceColor','k','LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1); % plot surface

    

    % Place electron on its orbit, including a spin vector arrow 

    elec = surf(r2*a-R,r2*b,r2*c,...

        'FaceColor',myBlue,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(elec,[0 0 1],elRot,[0,0,0]); 

    spinArrow = mArrow3([-R 0 -0.64],[-R 0 0.865],'color',myBlue,'stemWidth',0.04,...

        'tipWidth',0.08); 

    rotate(spinArrow,[0 0 1],elRot,[0,0,0]); 


    % Create a set of very thin torii to describe the magnetic-field lines

    for ii = 0:360/7:360

        for jj = 1:4

            R2=3; % outer radius of torus

            r=0.025; % inner tube radius

            th=linspace(0,2*pi,90); % e.g. 90 partitions along perimeter of the tube

            phi=linspace(0,2*pi,180); % e.g. 180 partitions along azimuth of torus

            % we convert our vectors phi and th to [n x n] matrices with meshgrid command:

            [Phi,Th]=meshgrid(phi,th);

            % now we generate n x n matrices for x,y,z according to eqn of torus

            x=(R2*1.4*jj/2+r.*cos(Th)).*cos(Phi);

            y=(R2*jj/2+r.*cos(Th)).*sin(Phi);

            z=r.*sin(Th);

            % Magnetic flux lines 

            magFL = surf(x,y-R+0.73-jj/1.1-jj^2/25,z,'FaceColor',myBlue,'LineStyle',...

                'none','FaceAlpha',magFalp,'FaceLighting','flat','DiffuseStrength',1); % plot surface

            rotate(magFL,[0 1 0],90,[0,0,0]);

            rotate(magFL,[0 0 1],ii,[0,0,0]);

            % Arrow head showing direction of magnetic flux lines 

            arrowHead = surf(X,Y-R+0.73+(cosd(elRot)*R2/2-1/1.1-jj/25)*jj,-sind(elRot)*R2*1.4*jj/2+Z*arrowLength-arrowLength/2,...

                'FaceColor',myBlue,'LineStyle',...

                'none','FaceLighting','flat','FaceAlpha',magFalp,'DiffuseStrength',1);

            rotate(arrowHead,[1 0 0],-elRot,[0,-R+0.73+(cosd(elRot)*R2/2-1/1.1-jj/25)*jj,-sind(elRot)*R2*1.4*jj/2]);            

            rotate(arrowHead,[0 0 1],ii,[0,0,0]);            

            arrowHeadDum = surf(X,Y-R+0.73+(cosd(180)*R2/2-1/1.1-jj/25)*jj,-sind(180)*R2*1.4*jj/2+Z*arrowLength-arrowLength/2,...

                'FaceColor',myBlue,'LineStyle',...

                'none','FaceLighting','flat','FaceAlpha',0,'DiffuseStrength',1);

            rotate(arrowHeadDum,[1 0 0],-180,[0,-R+0.73+(cosd(180)*R2/2-1/1.1-jj/25)*jj,-sind(180)*R2*1.4*jj/2]);            

            rotate(arrowHeadDum,[0 0 1],ii,[0,0,0]);            

        end

    end

    axis equal

    axis off

    axis tight

    light('Position',lp,'Style','infinite');

    set(gca, 'Projection','perspective');

    set(gca,'View',[139,50]); 

    

    

%__________________________________________________________________________ 

% Spin-orbit split levels j = l - s and j = l + s 

    subplot('position',pos2);

    newplot 

    hold off

    line([0 0.4],[0 0],'Color',[0.5 0.5 0.5],'LineWidth',5); % Unperturbed energy level

    hold on

    text(0.2,0.05,'l','FontSize',20,'HorizontalAlignment', 'center', ...

            'VerticalAlignment', 'bottom','color','k');

        

    % j = l - s energy level 

    rect = findall(gcf,'Type', 'Rectangle'); 

    delete(rect); 

    if (elRot >= 360)

        jlms = plot([0.6 1.0], [-1 -1], 'color',[0.5 0.5 0.5],'Linewidth', 5); 

        connect1 = plot([0.4 0.6], [0 -1], 'color',[0.75 0.75 0.75],'Linewidth', 2.5); 

    end

    if (elRot >= 360) && (elRot < 400)

        jlms.Color(4) = 1 - (400-elRot)/40; % Semitransparent line

        connect1.Color(4) = 1 - (400-elRot)/40; % Semitransparent line

        textCol = [1 1 1]*(400-elRot)/40;

        text(0.8,-0.95,'j = l - s','FontSize',20,'HorizontalAlignment', 'center', ...

            'VerticalAlignment', 'bottom','color',textCol);

        rectangle('Position',[0.775,-0.55,0.05,0.2],'FaceColor',[0 0 0 1-(400-elRot)/40],'EdgeColor','none')

        rectangle('Position',[0.775,-0.75,0.05,0.2],'FaceColor',[1 0 0 1-(400-elRot)/40],'EdgeColor','none')

        rectangle('Position',[0.85,-0.55,0.025,0.1],'FaceColor',[1 0 0 1-(400-elRot)/40],'EdgeColor','none')

        rectangle('Position',[0.85,-0.65,0.025,0.1],'FaceColor',[0 0 0 1-(400-elRot)/40],'EdgeColor','none')

    elseif(elRot >= 400) 

        text(0.8,-0.95,'j = l - s','FontSize',20,'HorizontalAlignment', 'center', ...

            'VerticalAlignment', 'bottom','color','k');

        rectangle('Position',[0.775,-0.55,0.05,0.2],'FaceColor',[0 0 0],'EdgeColor','none')

        rectangle('Position',[0.775,-0.75,0.05,0.2],'FaceColor',[1 0 0],'EdgeColor','none')

        rectangle('Position',[0.85,-0.55,0.025,0.1],'FaceColor',[1 0 0],'EdgeColor','none')

        rectangle('Position',[0.85,-0.65,0.025,0.1],'FaceColor',[0 0 0],'EdgeColor','none')

   end   

    

    % j = l + s energy level 

    if (elRot >= 720)

        jlps = line([0.6 1],[1 1],'Color',[0.5 0.5 0.5],'LineWidth',5);

        connect2 = plot([0.4 0.6], [0 1], 'color',[0.75 0.75 0.75],'Linewidth', 2.5); 

    end

    if (elRot >= 720) && (elRot < 760)

        jlps.Color(4) = 1 - (760-elRot)/40; % Semitransparent line

        connect2.Color(4) = 1 - (760-elRot)/40; % Semitransparent line

        textCol = [1 1 1]*(760-elRot)/40;

        text(0.8,1.05,'j = l + s','FontSize',20,'HorizontalAlignment', 'center', ...

            'VerticalAlignment', 'bottom','color',textCol);

        rectangle('Position',[0.77,1.45,0.05,0.2],'FaceColor',[0 0 0 1-(760-elRot)/40],'EdgeColor','none')

        rectangle('Position',[0.77,1.25,0.05,0.2],'FaceColor',[1 0 0 1-(760-elRot)/40],'EdgeColor','none')

        rectangle('Position',[0.85,1.35,0.025,0.1],'FaceColor',[1 0 0 1-(760-elRot)/40],'EdgeColor','none')

        rectangle('Position',[0.85,1.45,0.025,0.1],'FaceColor',[0 0 0 1-(760-elRot)/40],'EdgeColor','none')

    elseif(elRot >= 760) 

        text(0.8,1.05,'j = l + s','FontSize',20,'HorizontalAlignment', 'center', ...

            'VerticalAlignment', 'bottom','color','k');

        rectangle('Position',[0.77,1.45,0.05,0.2],'FaceColor',[0 0 0],'EdgeColor','none')

        rectangle('Position',[0.77,1.25,0.05,0.2],'FaceColor',[1 0 0],'EdgeColor','none')

        rectangle('Position',[0.85,1.35,0.025,0.1],'FaceColor',[1 0 0],'EdgeColor','none')

        rectangle('Position',[0.85,1.45,0.025,0.1],'FaceColor',[0 0 0],'EdgeColor','none')

    end

    axis off

    xlim([0 1]);

    ylim([-2 2]);

    

    frame = getframe(gcf);

    writeVideo(vid,frame);

end

% Output the movie as an mpg file

close(vid);