Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
Macros | Functions
nf_angularMomentumCoupling.cc File Reference
#include <stdlib.h>
#include <cmath>
#include "nf_specialFunctions.h"

Go to the source code of this file.

Macros

#define _USE_MATH_DEFINES
 

Functions

double nf_amc_log_factorial (int n)
 
double nf_amc_factorial (int n)
 
double nf_amc_wigner_3j (int j1, int j2, int j3, int j4, int j5, int j6)
 
double nf_amc_wigner_6j (int j1, int j2, int j3, int j4, int j5, int j6)
 
double nf_amc_wigner_9j (int j1, int j2, int j3, int j4, int j5, int j6, int j7, int j8, int j9)
 
double nf_amc_racah (int j1, int j2, int l2, int l1, int j3, int l3)
 
double nf_amc_clebsh_gordan (int j1, int j2, int m1, int m2, int j3)
 
double nf_amc_z_coefficient (int l1, int j1, int l2, int j2, int s, int ll)
 
double nf_amc_zbar_coefficient (int l1, int j1, int l2, int j2, int s, int ll)
 
double nf_amc_reduced_matrix_element (int lt, int st, int jt, int l0, int j0, int l1, int j1)
 

Macro Definition Documentation

◆ _USE_MATH_DEFINES

#define _USE_MATH_DEFINES

Definition at line 57 of file nf_angularMomentumCoupling.cc.

Function Documentation

◆ nf_amc_clebsh_gordan()

double nf_amc_clebsh_gordan ( int j1,
int j2,
int m1,
int m2,
int j3 )

Definition at line 289 of file nf_angularMomentumCoupling.cc.

289 {
290/*
291* Clebsh-Gordan coefficient
292* = <j1,j2,m1,m2|j3,m1+m2>
293* = (-)^(j1-j2+m1+m2) * std::sqrt(2*j3+1) * / j1 j2 j3 \
294* \ m1 m2 -m1-m2 /
295*
296* Note: Last value m3 is preset to m1+m2. Any other value will evaluate to 0.0.
297*/
298
299 int m3, x1, x2, x3, y1, y2, y3;
300 double cg = 0.0;
301
302 if ( j1 < 0 || j2 < 0 || j3 < 0) return( 0.0 );
303 if ( j1 + j2 + j3 > 2 * MAX_FACTORIAL ) return( INFINITY );
304
305 m3 = m1 + m2;
306
307 if ( ( x1 = ( j1 + m1 ) / 2 + 1 ) <= 0 ) return( 0.0 );
308 if ( ( x2 = ( j2 + m2 ) / 2 + 1 ) <= 0 ) return( 0.0 );
309 if ( ( x3 = ( j3 - m3 ) / 2 + 1 ) <= 0 ) return( 0.0 );
310
311 if ( ( y1 = x1 - m1 ) <= 0 ) return( 0.0 );
312 if ( ( y2 = x2 - m2 ) <= 0 ) return( 0.0 );
313 if ( ( y3 = x3 + m3 ) <= 0 ) return( 0.0 );
314
315 if ( j3 == 0 ){
316 if ( j1 == j2 ) cg = ( 1.0 / std::sqrt( (double)j1 + 1.0 ) * ( ( y1 % 2 == 0 ) ? -1:1 ) );
317 }
318 else if ( (j1 == 0 || j2 == 0 ) ){
319 if ( ( j1 + j2 ) == j3 ) cg = 1.0;
320 }
321 else {
322 if( m3 == 0 && std::abs( m1 ) <= 1 ){
323 if( m1 == 0 ) cg = cg1( x1, x2, x3 );
324 else cg = cg2( x1 + y1 - y2, x3 - 1, x1 + x2 - 2, x1 - y2, j1, j2, j3, m2 );
325 }
326 else if ( m2 == 0 && std::abs( m1 ) <=1 ){
327 cg = cg2( x1 - y2 + y3, x2 - 1, x1 + x3 - 2, x3 - y1, j1, j3, j3, m1 );
328 }
329 else if ( m1 == 0 && std::abs( m3 ) <= 1 ){
330 cg = cg2( x1, x1 - 1, x2 + x3 - 2, x2 - y3, j2, j3, j3, -m3 );
331 }
332 else cg = cg3( x1, x2, x3, y1, y2, y3 );
333 }
334
335 return( cg );
336}

Referenced by nf_amc_reduced_matrix_element(), nf_amc_wigner_3j(), nf_amc_z_coefficient(), and nf_amc_zbar_coefficient().

◆ nf_amc_factorial()

double nf_amc_factorial ( int n)

Definition at line 97 of file nf_angularMomentumCoupling.cc.

97 {
98/*
99* returns n! for pre-computed table. INFINITY is return if n is negative or too large.
100*/
101 return G4Exp( nf_amc_log_factorial( n ) );
102}
G4Exp
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition G4Exp.hh:180
nf_amc_log_factorial
double nf_amc_log_factorial(int n)
Definition nf_angularMomentumCoupling.cc:86

◆ nf_amc_log_factorial()

double nf_amc_log_factorial ( int n)

Definition at line 86 of file nf_angularMomentumCoupling.cc.

86 {
87/*
88* returns ln( n! ).
89*/
90 if( n > MAX_FACTORIAL ) return( INFINITY );
91 if( n < 0 ) return( INFINITY );
92 return nf_amc_log_fact[n];
93}

Referenced by nf_amc_factorial().

◆ nf_amc_racah()

double nf_amc_racah ( int j1,
int j2,
int l2,
int l1,
int j3,
int l3 )

Definition at line 254 of file nf_angularMomentumCoupling.cc.

254 {
255/*
256* Racah coefficient definition in Edmonds (AR Edmonds, "Angular Momentum in Quantum Mechanics", Princeton (1980) is
257* W(j1, j2, l2, l1 ; j3, l3) = (-1)^(j1+j2+l1+l2) * { j1 j2 j3 }
258* { l1 l2 l3 }
259* The call signature of W(...) appears jumbled, but hey, that's the convention.
260*
261* This convention is exactly that used by Blatt-Biedenharn (Rev. Mod. Phys. 24, 258 (1952)) too
262*/
263
264 double sig;
265
266 sig = ( ( ( j1 + j2 + l1 + l2 ) % 4 == 0 ) ? 1.0 : -1.0 );
267 return sig * nf_amc_wigner_6j( j1, j2, j3, l1, l2, l3 );
268}
nf_amc_wigner_6j
double nf_amc_wigner_6j(int j1, int j2, int j3, int j4, int j5, int j6)
Definition nf_angularMomentumCoupling.cc:122

Referenced by nf_amc_wigner_9j(), nf_amc_z_coefficient(), and nf_amc_zbar_coefficient().

◆ nf_amc_reduced_matrix_element()

double nf_amc_reduced_matrix_element ( int lt,
int st,
int jt,
int l0,
int j0,
int l1,
int j1 )

Definition at line 474 of file nf_angularMomentumCoupling.cc.

474 {
475/*
476* Reduced Matrix Element for Tensor Operator
477* = < l1j1 || T(YL,sigma_S)J || l0j0 >
478*
479* M.B.Johnson, L.W.Owen, G.R.Satchler
480* Phys. Rev. 142, 748 (1966)
481* Note: definition differs from JOS by the factor sqrt(2j1+1)
482*/
483 int llt;
484 double x1, x2, x3, reduced_mat, clebsh_gordan;
485
486 if ( parity( lt ) != parity( l0 ) * parity( l1 ) ) return( 0.0 );
487 if ( std::abs( l0 - l1 ) > lt || ( l0 + l1 ) < lt ) return( 0.0 );
488 if ( std::abs( ( j0 - j1 ) / 2 ) > jt || ( ( j0 + j1 ) / 2 ) < jt ) return( 0.0 );
489
490 llt = 2 * lt;
491 jt *= 2;
492 st *= 2;
493
494 if( ( clebsh_gordan = nf_amc_clebsh_gordan( j1, j0, 1, -1, jt ) ) == INFINITY ) return( INFINITY );
495
496 reduced_mat = 1.0 / std::sqrt( 4 * M_PI ) * clebsh_gordan / std::sqrt( jt + 1.0 ) /* BRB jt + 1 <= 0? */
497 * std::sqrt( ( j0 + 1.0 ) * ( j1 + 1.0 ) * ( llt + 1.0 ) )
498 * parity( ( j1 - j0 ) / 2 ) * parity( ( -l0 + l1 + lt ) / 2 ) * parity( ( j0 - 1 ) / 2 );
499
500 if( st == 2 ){
501 x1 = ( l0 - j0 / 2.0 ) * ( j0 + 1.0 );
502 x2 = ( l1 - j1 / 2.0 ) * ( j1 + 1.0 );
503 if ( jt == llt ){
504 x3 = ( lt == 0 ) ? 0 : ( x1 - x2 ) / std::sqrt( lt * ( lt + 1.0 ) );
505 }
506 else if ( jt == ( llt - st ) ){
507 x3 = ( lt == 0 ) ? 0 : -( lt + x1 + x2 ) / std::sqrt( lt * ( 2.0 * lt + 1.0 ) );
508 }
509 else if ( jt == ( llt + st ) ){
510 x3 = ( lt + 1 - x1 - x2 ) / std::sqrt( ( 2.0 * lt + 1.0 ) * ( lt + 1.0 ) );
511 }
512 else{
513 x3 = 1.0;
514 }
515 }
516 else x3 = 1.0;
517 reduced_mat *= x3;
518
519 return( reduced_mat );
520}
M_PI
#define M_PI
Definition SbMath.h:33
nf_amc_clebsh_gordan
double nf_amc_clebsh_gordan(int j1, int j2, int m1, int m2, int j3)
Definition nf_angularMomentumCoupling.cc:289

◆ nf_amc_wigner_3j()

double nf_amc_wigner_3j ( int j1,
int j2,
int j3,
int j4,
int j5,
int j6 )

Definition at line 106 of file nf_angularMomentumCoupling.cc.

106 {
107/*
108* Wigner's 3J symbol (similar to Clebsh-Gordan)
109* = / j1 j2 j3 \
110* \ j4 j5 j6 /
111*/
112 double cg;
113
114 if( ( j4 + j5 + j6 ) != 0 ) return( 0.0 );
115 if( ( cg = nf_amc_clebsh_gordan( j1, j2, j4, j5, j3 ) ) == 0.0 ) return ( 0.0 );
116 if( cg == INFINITY ) return( cg );
117 return( ( ( ( j1 - j2 - j6 ) % 4 == 0 ) ? 1.0 : -1.0 ) * cg / std::sqrt( j3 + 1.0 ) ); /* BRB j3 + 1 <= 0? */
118}

◆ nf_amc_wigner_6j()

double nf_amc_wigner_6j ( int j1,
int j2,
int j3,
int j4,
int j5,
int j6 )

Definition at line 122 of file nf_angularMomentumCoupling.cc.

122 {
123/*
124* Wigner's 6J symbol (similar to Racah)
125* = { j1 j2 j3 }
126* { j4 j5 j6 }
127*/
128 int i, x[6];
129
130 x[0] = j1; x[1] = j2; x[2] = j3; x[3] = j4; x[4] = j5; x[5] = j6;
131 for( i = 0; i < 6; i++ ) if ( x[i] == 0 ) return( w6j0( i, x ) );
132
133 return( w6j1( x ) );
134}

Referenced by nf_amc_racah().

◆ nf_amc_wigner_9j()

double nf_amc_wigner_9j ( int j1,
int j2,
int j3,
int j4,
int j5,
int j6,
int j7,
int j8,
int j9 )

Definition at line 227 of file nf_angularMomentumCoupling.cc.

227 {
228/*
229* Wigner's 9J symbol
230* / j1 j2 j3 \
231* = | j4 j5 j6 |
232* \ j7 j8 j9 /
233*
234*/
235 int i, i0, i1;
236 double rac;
237
238 i0 = max3( std::abs( j1 - j9 ), std::abs( j2 - j6 ), std::abs( j4 - j8 ) );
239 i1 = min3( ( j1 + j9 ), ( j2 + j6 ), ( j4 + j8 ) );
240
241 rac = 0.0;
242 for ( i = i0; i <= i1; i += 2 ){
243 rac += nf_amc_racah( j1, j4, j9, j8, j7, i )
244 * nf_amc_racah( j2, j5, i, j4, j8, j6 )
245 * nf_amc_racah( j9, i, j3, j2, j1, j6 ) * ( i + 1 );
246 if( rac == INFINITY ) return( INFINITY );
247 }
248
249 return( ( ( (int)( ( j1 + j3 + j5 + j8 ) / 2 + j2 + j4 + j9 ) % 4 == 0 ) ? 1.0 : -1.0 ) * rac );
250}
nf_amc_racah
double nf_amc_racah(int j1, int j2, int l2, int l1, int j3, int l3)
Definition nf_angularMomentumCoupling.cc:254

◆ nf_amc_z_coefficient()

double nf_amc_z_coefficient ( int l1,
int j1,
int l2,
int j2,
int s,
int ll )

Definition at line 438 of file nf_angularMomentumCoupling.cc.

438 {
439/*
440* Biedenharn's Z-coefficient coefficient
441* = Z(l1 j1 l2 j2 | S L )
442*/
443 double z, clebsh_gordan = nf_amc_clebsh_gordan( l1, l2, 0, 0, ll ), racah = nf_amc_racah( l1, j1, l2, j2, s, ll );
444
445 if( ( clebsh_gordan == INFINITY ) || ( racah == INFINITY ) ) return( INFINITY );
446 z = ( ( ( -l1 + l2 + ll ) % 8 == 0 ) ? 1.0 : -1.0 )
447 * std::sqrt( l1 + 1.0 ) * std::sqrt( l2 + 1.0 ) * std::sqrt( j1 + 1.0 ) * std::sqrt( j2 + 1.0 ) * clebsh_gordan * racah;
448
449 return( z );
450}

◆ nf_amc_zbar_coefficient()

double nf_amc_zbar_coefficient ( int l1,
int j1,
int l2,
int j2,
int s,
int ll )

Definition at line 454 of file nf_angularMomentumCoupling.cc.

454 {
455/*
456* Lane & Thomas's Zbar-coefficient coefficient
457* = Zbar(l1 j1 l2 j2 | S L )
458* = (-i)^( -l1 + l2 + ll ) * Z(l1 j1 l2 j2 | S L )
459*
460* Lane & Thomas Rev. Mod. Phys. 30, 257-353 (1958).
461* Note, Lane & Thomas define this because they did not like the different phase convention in Blatt & Biedenharn's Z coefficient. They changed it to get better time-reversal behavior.
462* Froehner uses Lane & Thomas convention as does T. Kawano.
463*/
464 double zbar, clebsh_gordan = nf_amc_clebsh_gordan( l1, l2, 0, 0, ll ), racah = nf_amc_racah( l1, j1, l2, j2, s, ll );
465
466 if( ( clebsh_gordan == INFINITY ) || ( racah == INFINITY ) ) return( INFINITY );
467 zbar = std::sqrt( l1 + 1.0 ) * std::sqrt( l2 + 1.0 ) * std::sqrt( j1 + 1.0 ) * std::sqrt( j2 + 1.0 ) * clebsh_gordan * racah;
468
469 return( zbar );
470}