This is NSSS version 0.1.0, a N-Spin System Simulator (NSSS) for NMR that is fully based on the Cartesian Operator Product Formalism.
This software was initially written (in 1994) because Ray Freeman asked me if I could check his cartesian operator product calculation for the DEPT sequence.
My second intention was to write some code for a virtual spectrometer, something I never did.
NSSS is a C API. The user has to write a C program to use it. Demo executable programs named dept and abc1d can be built from the distribution file. They provide examples of NSSS use.
typedef unsigned int usint ;
typedef struct Term
{
the address of the term to which this term is linked
struct Term *next ;
this flag indicated whether this term corresponds to a state (value: STAT)
or to an operator (value: HAMI). Also used to store a pointer
in base and commuter terms. see subspace().
int flag ;
loc1 and loc2 are two multi-purpose integer values.
int loc1, loc2 ;
the term axis set
usint axis ;
the factor that multiplies the operator part of the term
double fact ;
}
A "term" structure stores something like "0.8 IzJxSz".
The axis set "IzJxSz" is coded in a single unsigned integer (usint),
that must be represented by 32 bits. This is not directly
portable to machines for which an "int" has less then 32 bits.
An axis in coded by 2 bits, for 4 values " " (nothing), "x", "y" and "z".
This limits the number of spins to 16.
MAXSPIN is set to 15 to avoid potential problems.
NSSS was never tester with more than 4 spins.
typedef term *pterm ;
typedef int (*pfunc)(pterm) ;
typedef struct s3FCOMPLEX
{
double r ;
double i ;
} s3fcomplex ;
A s3fcomplex is simply a complex number.
All operator algebra in NSSS is performed with real (double) numbers.
Complexe numbers are needed for FID and phase modulated shaped pulses
manipulation.
#define HAMI 1
#define STAT 0
define two types of terms, HAMI for operators and STAT for spin states.
#define SHORT 0
#define LONG 1
#define PRE2C 1.0e-9
define two ways of printing a sum of terms.
The LONG version prints all present terms.
The SHORT one prints only term for which
abs(factor) is bigger than PRE2C.
void spins(char *s)
declares s as spin system. Argument looks like "IJS".
Spin identifiers ('I', 'J', 'S') are capital letters.
spins() sets letter_to_num so that letter_to_num[9] is 1 because 'J'
is the 10th letter of the alphabet (index 9) and the 2nd (index 1)
letter in the "IJS" string.
Then, num_to_letter[1] is 'J'.
spins also sets nspins, the number of spins in the system.
pterm create_term(pterm anchor)
creates a new term and links it to its anchor. anchor must initially
be NULL, else your program will crash.
usint desc_to_axis(char *desc)
builds an axis set representation
from desc, its string description,
something that looks like "IxJySz".
void fill_term(pterm p, char *desc, double f, int typ)
sets fields in an already existing term pointed by p.
desc is the axis set descriptor, f the factor, and typ
is HAMI or STAT.
pterm create_state(pterm anchor, char *desc, double f)
creates a state term of descriptor desc, factor f, and anchors
it at anchor. Returns its address as new anchor for an
eventualy following term.
pterm create_hami(pterm anchor, char *desc, double f)
creates a new Hamiltonian term of descriptor desc, factor f in rad/s,
and anchors it at anchor. Returns its address as new anchor for an
eventually additionnal term.
pterm duplicate_state(pterm anchor)
returns an anchor to a newly created sum of STAT terms
that carries the same information as anchor.
void modif_term(pterm anchor, char *desc, double nf)
changes to nf the factor that is associated to a axis description desc
within a term sum that is pointed by anchor.
void write_term(pterm p)
writes a single term that is pointed by p.
void write_term_with_way(pterm p, int way)
writes a single term, always if way is LONG and when abs(fact) is
bigger than the PRE2C threshold if way is SHORT.
void write_terms(pterm anchor, int way)
writes a sum of terms according to way (SHORT or LONG).
void write_terms_if(pterm anchor, int way, pfunc f)
writes a sum of terms according to way (SHORT or LONG).
The terms for which function f() returns 0 are not written.
int quanta(pterm p)
returns the highest m_s(total) quantum number value of the
single term that is pointed by p.
void detected(char *s)
defines by character string s the spin set that will be detected.
int detectable(pterm p)
returns the detectability of a term. 0 means not detectable.
1 stands for an X state and 2 for a Y state.
s3fcomplex detect(pterm anchor)
computes the signal that arises from the detection of a sum of terms
that is pointed by anchor. It uses previous argument to detected()
to select which nuclei are detected.
pterm extract(pterm anchor, pfunc f)
in-place removes terms in a term sum pointed by anchor
that do not satisfy the f() function.
int is_single_quantum(pterm p)
returns the single_quantumness of a term.
pterm single_quantum(pterm anchor)
keeps in a term sum pointed by anchor
only the single quantum ones.
pterm sum_term(pterm anchor, usint axes, double f)
adds a single new term defined by axes and factor to the sum
that is pointed by anchor.
pterm sum_terms(pterm buf, pterm anchor, double f)
accumulates terms from anchor in buf, returns the address of the sum.
void comm(pterm p, pterm h, usint *na, int *cs)
calculates the axis set and sign of the {h, p} operator commuter.
na points to the new axis, and cs to the commuter sign.
*cs and *na are zero if h and p commute
pterm action_hami(pterm anchor, pterm h, double t)
computes the result of the application of a single hamiltonian term h
on a sum of terms starting at anchor.
pterm action(pterm anchor, pterm h, double t)
implements product of cartesian operator algebra. Hamiltonian h is applied
to state anchor during time t. action returns the final state of the system.
Each term in h must commute with all others, otherwise unpredictable (wrong)
results will be produced. This is not checked for.
pterm simplify(pterm anchor)
removes terms with a zero factor from a term sum.
int length(pterm anchor)
returns the number of terms in a term sum
pterm look_axis(pterm anchor, usint a)
returns the address of the first term in sum anchor that has a as axis.
pterm look_desc(pterm anchor, char *desc)
returns the address of the first term in sum anchor that has desc
as axis set descriptor
int desc_to_index(pterm anchor, char *desc)
returns the index within the base pointed by anchor,
of a term with desc as axis set descriptor
void subspace(pterm anchor, pterm h, pterm *base, int *size)
returns a pointer to a base and its size(numner of elements).
The base contains at least all the terms pointed by anchor and is stable
by action of operator sum h. This means that the action of h on
any linear combination of base terms is a linear combination of
base terms. This also means that the action of h on anchor
leads to a state that is a linear combination of the operators in base.
Subspace() must be called only once for a given h and a
given anchor. Subspace() does not compute anything. It only sets
pointers so that bch() can transform an initial state with the
same axis sets than the one pointed by anchor by action of h.
See the code of write_fid() for an example of how this can be used.
pterm to_be_traced(pterm anchor, char *desc)
builds a list of operators that will be used
for trace (printout) generation. See init_trace().
The problem is to print the factors of a term sum
in the desired order, independently of the order
in which they are stored. The desired order is defined
by a sum of terms that is build by successive calls
to to_be_traced().
int *init_trace(pterm anchor, pterm base)
returns a newly created trace vector from a list of terms to be traced,
pointed by anchor, so that the factors of the terms in base
(see subspace()) will correctly be printed by trace_vec(), once
tranferred from a base to a vector by state_to_vec().
Any sum of terms may be considered as a base by init_trace()
after it has been processed by itemize().
void trace_vec(int *tv, double *v)
prints the elements in v (generated by state_to_vec()) according to
indexes in tv that were computed by init_trace().
void itemize(pterm anchor)
initializes loc2 in a sum of terms so that desc_to_index(),
init_trace(), state_to_vec can work.
int size_term(pterm anchor)
looks very much like length() -- should probably be removed.
Kept for compatibility with old programs
double *state_to_vec(int size, pterm base, pterm anchor)
returns a pointer to a new vector that describes the anchor state
in the base base of length size.
void vec_to_base(int size, double *v, pterm base)
stores the size vector elements in v as factors in base
void bch(int size, double *ini, double *res, pterm base, pterm h, double t)
transforms the ini initial state vector of length size,
expressed in base base, when submitted to hamiltonian h during t seconds.
The resulting vector is stored at localtion res.
Of course, subspace() must have previously called to define base from h
and an initial state that looks like the one that was transformed
by state_to_vec() to give ini. Looks like means same operators but
possibly different factors
bch computes (exp(-iH''t))(ini) in the base using Taylor expansion.
H'' (H "double hat") is the superoperator that is is associated
to operator H' (H "single hat") that is defined by H''(A') = {H', A'}.
The Baker-Campbell-Hausdorff (or BCH) formula is simply the developped
expression of the Taylor expansion.
double **comm_mat(int size, pterm base)
computes the matrix of the commutation superoperator
of the hamiltonian that what used by subspace()
to initialize base. The length of base is size.
Comm_mat returns the result as a newly created matrix
made of size lines and size columns.
void matrot(int size, double **m, pterm base, pterm h, double t)
computes m in base base, the matrix of the superoperator that corresponds to
the action of hamiltonian h during time t. base must have been previously
initialized using h as argument to subspace().
void init_fid(char *nuc, pterm hami, int np, double sw, double T2)
initializes an FID computation for a nuclei set nuc (calls detected(nuc)),
a free evolution hamiltonian h, a number of complex points np,
a spectral width sw, and an apparent tranverse relaxation time T2.
A negative T2 is considered as an infinite T2.
void print_fid(pterm state, float *buf)
does not print anything but computes the FID that is generated
by the evolution of the initial state state, according to
parameters given to init_fid. The FID is stored at location buf,
whose size is at least 2*np. The data are stored as
real0, imag0, real1, imag1, real2, ... for data points 0, 1, 2...
The code of this function shows how most of the Simu3 API functions
work together.
void copy_vector(int size, double *org, double *dest)
copies a vector of floats from org to dest
void display_mat(int size, double **m)
displays the terms of the square matrix m, that has size lines and columns.
void unite(int size, double **m)
sets an existing square matrix of size size to the unit matrix.
void matmult(int size, double **m, double *v1, double *v2)
computes v2 as the product of matrix v by vector v1, all of
them being of size size. This function should be renamed matvec()
to be consistent with libsimu1 and libsimu1b libraries.
void matmat(int size, double **m1, double **m2, double **m3)
computes matrix m3 as the product of matrices m2 and m3,
all of them being of size size.
s3fcomplex s3Cadd(s3fcomplex a, s3fcomplex b)
returns of sum of two complexe numbers.
s3fcomplex s3Csub(s3fcomplex a, s3fcomplex b)
returns of difference of two complexe numbers.
s3fcomplex s3CxD(s3fcomplex a, double b)
returns of product of a complexe number a by a double b.
s3fcomplex s3Complex(double re, double im)
makes a complexe number from two doubles.
double **s3dmatrix(int nr, int nc)
returns a matrix of doubles with nr rows and nc columns.
void s3free_dmatrix(double **m)
frees the storage space of a matrix of doubles.
float **s3fmatrix(int nr, int nc)
returns a matrix of floats with nr rows and nc columns.
void s3free_fmatrix(float **m)
frees the storage space of a matrix of floats.
double *s3dvector(int n)
returns a vector of doubles with n elements.
void s3free_dvector(double *v)
frees the storage space of a vector of doubles.
float *s3fvector(int n)
returns a vector of floats with n elements.
void s3free_fvector(float *v)
frees the storage space of a vector of floats.
int *s3ivector(int n)
returns a vector of integers with n elements.
s3fcomplex *s3cvector(int n)
returns a vector of complexes with n elements.
void entering(char *s)
keeps track of the name of the last called function, see s3erreur().
A true stack of function calls would be more useful, may be.
void s3erreur(char *s)
prints error message s and the name of the last
function that called entering()
Jean-Marc Nuzillard
jm.nuzillard@univ-reims.fr
April 11^{th} 2006