src/HSLib/SpinOpCmp.cc
/* The spin operators returned here reside in the composite Hilbert space
   of the spin system given in the function call even though they operate
   only on a single spin!  For example, Ix for a single proton would be
   represented by a 2x2 array whereas that same protons Ix when part of
   a three proton system would by 8x8 where 8 is the system Hilbert space.

		Input		sys	: A spin system
				spin	: A spin index [0, ns)
		Output		SOp	: A spin operator in the full
					  Hilbert space of sys               */
The spin operators returned here reside in the composite Hilbert space
of the spin system given in the function call even though they operate
only on a single spin.  For example, Ix for a single proton would be
represented by a 2x2 array whereas that same protons Ix when part of
a three proton system would by 8x8 where 8 is the system Hilbert space.

Input           sys     : A spin system
spin    : A spin index [0, ns)
Output          SOp     : A spin operator in the full
Hilbert space of sys
The spin operators returned here reside in the composite Hilbert space of the spin system given in the function call even though they operate only on a single spin . For example Ix for a single proton would be represented by a 2x2 array whereas that same protons Ix when part of a three proton system would by 8x8 where 8 is the system Hilbert space . Input sys : A spin system spin : A spin index [ 0 ns ) Output SOp : A spin operator in the full Hilbert space of sys
spin_op Iu(const spin_sys& sys, int spin, int type)
  {
  int hs, ns=sys.spins();		
  matrix* spr = new matrix[ns];		
  for(int i=0; i<ns; i++)		
    {
    hs = sys.HS(i);			
    if(spin == i)			
      {					
      switch(type)
        {
        default:
        case 0: spr[i]=Ie(hs); break;	
        case 1: spr[i]=Ix(hs); break;	
        case 2: spr[i]=Iy(hs); break;	
        case 3: spr[i]=Iz(hs); break;	
        case 4: spr[i]=Ip(hs); break;	
        case 5: spr[i]=Im(hs); break;	
        case 6:				
          if(hs != 2) 			
            {				
            SOpCmperror(20, 1);		
            SOpCmpfatal(22);		
            }
          spr[i]=matrix(hs,hs,0,d_matrix_type);
          (spr[i]).put(1.0, 0,0);
          break;
        case 7:				
          if(hs != 2) 			
            {				
            SOpCmperror(21, 1);		
            SOpCmpfatal(22);		
            }
          spr[i]=matrix(hs,hs,0,d_matrix_type);
          (spr[i]).put(1.0, 1,1);
          break;
        }
      }
    else spr[i] = Ie(hs);		
    }
  spin_op SOp(ns, spr);
  if(spr) delete [] spr;
  return SOp;
  }
spin_op Iu ( const spin_sys & sys , int spin , int type ) { int hs , ns = sys . spins ( ) ; matrix * spr = new matrix [ ns ] ; for ( int i = 0 ; i < ns ; i ++ ) { hs = sys . HS ( i ) ; if ( spin == i ) { switch ( type ) { default : case 0 : spr [ i ] = Ie ( hs ) ; break ; case 1 : spr [ i ] = Ix ( hs ) ; break ; case 2 : spr [ i ] = Iy ( hs ) ; break ; case 3 : spr [ i ] = Iz ( hs ) ; break ; case 4 : spr [ i ] = Ip ( hs ) ; break ; case 5 : spr [ i ] = Im ( hs ) ; break ; case 6 : if ( hs != 2 ) { SOpCmperror ( 20 , 1 ) ; SOpCmpfatal ( 22 ) ; } spr [ i ] = matrix ( hs , hs , 0 , d_matrix_type ) ; ( spr [ i ] ) . put ( 1.0 , 0 , 0 ) ; break ; case 7 : if ( hs != 2 ) { SOpCmperror ( 21 , 1 ) ; SOpCmpfatal ( 22 ) ; } spr [ i ] = matrix ( hs , hs , 0 , d_matrix_type ) ; ( spr [ i ] ) . put ( 1.0 , 1 , 1 ) ; break ; } } else spr [ i ] = Ie ( hs ) ; } spin_op SOp ( ns , spr ) ; if ( spr ) delete [ ] spr ; return SOp ; }
The spin operators returned here reside in the composite Hilbert space
of the spin system given in the function call even though they operate
only on a single spin!
The spin operators returned here reside in the composite Hilbert space of the spin system given in the function call even though they operate only on a single spin!
// Number of system spins // Array of sub-space Ie's // Loop over all spins //	Spin i Hilbert space //	If this is the active spin //	we must make proper sub-matirx //	For Ie operator //	For Ix operator //	For Iy operator //	For Iz operator //	For I+ operator //	For I- operator //	For Ialpha operator // 	Error: Ia on I>1/2 spin //	Should use Ipol instead //      Cannot have Ia on I>1/2 //	Died on polarization Op //	For Ibeta operator // 	Error: Ib on I>1/2 spin //	Should use Ipol instead //      Cannot have Ib on I>1/2 //	Died on polarization Op // If not active spin, just Ie