Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
nf_Legendre.h File Reference
#include <nf_utilities.h>
#include <ptwXY.h>
Include dependency graph for nf_Legendre.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Classes

struct  nf_Legendre_s
 

Macros

#define nf_Legendre_minMaxOrder   4
 
#define nf_Legendre_maxMaxOrder   64
 
#define nf_Legendre_sizeIncrement   8
 

Typedefs

typedef struct nf_Legendre_s nf_Legendre
 
typedef nfu_status(* nf_Legendre_GaussianQuadrature_callback )(double x, double *y, void *argList)
 

Functions

nf_Legendrenf_Legendre_new (int initialSize, int maxOrder, double *Cls, nfu_status *status)
 
nfu_status nf_Legendre_setup (nf_Legendre *nfL, int initialSize, int maxOrder)
 
nfu_status nf_Legendre_release (nf_Legendre *nfL)
 
nf_Legendrenf_Legendre_free (nf_Legendre *nfL)
 
nf_Legendrenf_Legendre_clone (nf_Legendre *nfL, nfu_status *status)
 
nfu_status nf_Legendre_reallocateCls (nf_Legendre *Legendre, int size, int forceSmallerResize)
 
int nf_Legendre_maxOrder (nf_Legendre *Legendre)
 
int nf_Legendre_allocated (nf_Legendre *Legendre)
 
double nf_Legendre_getCl (nf_Legendre *Legendre, int l, nfu_status *status)
 
nfu_status nf_Legendre_setCl (nf_Legendre *Legendre, int l, double Cl)
 
nfu_status nf_Legendre_normalize (nf_Legendre *Legendre)
 
double nf_Legendre_evauluateAtMu (nf_Legendre *nfL, double mu, nfu_status *status)
 
double nf_Legendre_PofL_atMu (int l, double mu)
 
ptwXYPointsnf_Legendre_to_ptwXY (nf_Legendre *nfL, double accuracy, int biSectionMax, int checkForRoots, nfu_status *status)
 
nf_Legendrenf_Legendre_from_ptwXY (ptwXYPoints *ptwXY, int maxOrder, nfu_status *status)
 
nfu_status nf_Legendre_GaussianQuadrature (int degree, double x1, double x2, nf_Legendre_GaussianQuadrature_callback func, void *argList, double *integral)
 

Macro Definition Documentation

#define nf_Legendre_maxMaxOrder   64

Definition at line 18 of file nf_Legendre.h.

#define nf_Legendre_minMaxOrder   4

Definition at line 17 of file nf_Legendre.h.

#define nf_Legendre_sizeIncrement   8

Definition at line 19 of file nf_Legendre.h.

Typedef Documentation

typedef struct nf_Legendre_s nf_Legendre

Definition at line 21 of file nf_Legendre.h.

typedef nfu_status(* nf_Legendre_GaussianQuadrature_callback)(double x, double *y, void *argList)

Definition at line 29 of file nf_Legendre.h.

Function Documentation

int nf_Legendre_allocated ( nf_Legendre Legendre)

Definition at line 112 of file nf_Legendre.cc.

112  {
113 
114  return( Legendre->allocated );
115 }
nf_Legendre* nf_Legendre_clone ( nf_Legendre nfL,
nfu_status status 
)

Definition at line 70 of file nf_Legendre.cc.

70  {
71 
72  return( nf_Legendre_new( 0, nfL->maxOrder, nfL->Cls, status ) );
73 }
nf_Legendre * nf_Legendre_new(int initialSize, int maxOrder, double *Cls, nfu_status *status)
Definition: nf_Legendre.cc:23
double * Cls
Definition: nf_Legendre.h:26

Here is the call graph for this function:

double nf_Legendre_evauluateAtMu ( nf_Legendre nfL,
double  mu,
nfu_status status 
)

Definition at line 160 of file nf_Legendre.cc.

160  {
161 
162  int l;
163  double P = 0.;
164 
165  *status = nfu_XOutsideDomain;
166  if( ( mu >= -1. ) && ( mu <= 1. ) ) {
167  *status = nfu_Okay;
168  for( l = 0; l <= Legendre->maxOrder; l++ ) P += ( l + 0.5 ) * Legendre->Cls[l] * nf_Legendre_PofL_atMu( l, mu );
169  }
170  return( P );
171 }
double nf_Legendre_PofL_atMu(int l, double mu)
Definition: nf_Legendre.cc:175
static double P[]

Here is the call graph for this function:

Here is the caller graph for this function:

nf_Legendre* nf_Legendre_free ( nf_Legendre nfL)

Definition at line 61 of file nf_Legendre.cc.

61  {
62 
63  nf_Legendre_release( Legendre );
64  nfu_free( Legendre );
65  return( NULL );
66 }
nfu_status nf_Legendre_release(nf_Legendre *nfL)
Definition: nf_Legendre.cc:52
void * nfu_free(void *p)

Here is the call graph for this function:

Here is the caller graph for this function:

nf_Legendre* nf_Legendre_from_ptwXY ( ptwXYPoints ptwXY,
int  maxOrder,
nfu_status status 
)

Definition at line 250 of file nf_Legendre.cc.

250  {
251 
252  int l, i, n = (int) ptwXY_length( ptwXY );
253  nf_Legendre *Legendre;
254  double mu1, mu2, f1, f2, Cl, Cls[1] = { 0 }, integral;
255  struct nf_Legendre_from_ptwXY_callback_s argList;
256 
257  if( ( *status = ptwXY_getStatus( ptwXY ) ) != nfu_Okay ) return( NULL );
258 
259  ptwXY_getXYPairAtIndex( ptwXY, 0, &mu1, &f1 );
260  if( mu1 < -1 ) {
261  *status = nfu_XOutsideDomain;
262  return( NULL );
263  }
264 
265  ptwXY_getXYPairAtIndex( ptwXY, 0, &mu2, &f2 );
266  if( mu2 > 1 ) {
267  *status = nfu_XOutsideDomain;
268  return( NULL );
269  }
270 
271  if( ( Legendre = nf_Legendre_new( maxOrder + 1, -1, Cls, status ) ) == NULL ) return( NULL );
272 
273  if( maxOrder > nf_Legendre_maxMaxOrder ) maxOrder = nf_Legendre_maxMaxOrder;
274  for( l = 0; l <= maxOrder; l++ ) {
275  ptwXY_getXYPairAtIndex( ptwXY, 0, &mu1, &f1 );
276  argList.l = l;
277  for( i = 1, Cl = 0; i < n; i++ ) {
278  ptwXY_getXYPairAtIndex( ptwXY, i, &mu2, &f2 );
279  argList.mu1 = mu1;
280  argList.f1 = f1;
281  argList.mu2 = mu2;
282  argList.f2 = f2;
283  if( ( *status = nf_Legendre_GaussianQuadrature( l + 1, mu1, mu2, nf_Legendre_from_ptwXY_callback, (void *) &argList, &integral ) ) != nfu_Okay )
284  goto err;
285  Cl += integral;
286  mu1 = mu2;
287  f1 = f2;
288  }
289  if( ( *status = nf_Legendre_setCl( Legendre, l, Cl ) ) != nfu_Okay ) goto err;
290  }
291  return( Legendre );
292 
293 err:
294  nf_Legendre_free( Legendre );
295  return( NULL );
296 }
nf_Legendre * nf_Legendre_new(int initialSize, int maxOrder, double *Cls, nfu_status *status)
Definition: nf_Legendre.cc:23
nfu_status nf_Legendre_setCl(nf_Legendre *Legendre, int l, double Cl)
Definition: nf_Legendre.cc:131
int64_t ptwXY_length(ptwXYPoints *ptwXY)
Definition: ptwXY_core.cc:583
#define nf_Legendre_maxMaxOrder
Definition: nf_Legendre.h:18
nfu_status ptwXY_getStatus(ptwXYPoints *ptwXY)
Definition: ptwXY_core.cc:351
typedef int(XMLCALL *XML_NotStandaloneHandler)(void *userData)
const G4int n
nf_Legendre * nf_Legendre_free(nf_Legendre *nfL)
Definition: nf_Legendre.cc:61
static nfu_status nf_Legendre_from_ptwXY_callback(double mu, double *f, void *argList)
Definition: nf_Legendre.cc:300
nfu_status nf_Legendre_GaussianQuadrature(int degree, double x1, double x2, nf_Legendre_GaussianQuadrature_callback func, void *argList, double *integral)
nfu_status ptwXY_getXYPairAtIndex(ptwXYPoints *ptwXY, int64_t index, double *x, double *y)
Definition: ptwXY_core.cc:698

Here is the call graph for this function:

nfu_status nf_Legendre_GaussianQuadrature ( int  degree,
double  x1,
double  x2,
nf_Legendre_GaussianQuadrature_callback  func,
void argList,
double *  integral 
)

Definition at line 63 of file nf_Legendre_GaussianQuadrature.cc.

63  {
64 
65  int i, n;
66  double x, mu, sum, *weights, *xis;
67  nfu_status status = nfu_Okay;
68 
69  *integral = 0;
70  if( degree < 2 ) {
71  status = func( 0.5 * ( x1 + x2 ), integral, argList );
72  *integral *= 2.; }
73  else if( degree < 4 ) {
74  x = 0.5 * ( -sqrt_inv3 * ( x2 - x1 ) + x1 + x2 );
75  if( ( status = func( x, integral, argList ) ) == nfu_Okay ) {
76  x = 0.5 * ( sqrt_inv3 * ( x2 - x1 ) + x1 + x2 );
77  status = func( x, &sum, argList );
78  *integral += sum;
79  } }
80  else {
81  for( i = 0; i < nSets - 1; i++ ) {
82  if( GaussianQuadrature_degrees[i].n > ( degree + 1 ) / 2 ) break;
83  }
84  n = ( GaussianQuadrature_degrees[i].n + 1 ) / 2;
87  for( i = 0; i < n; i++ ) {
88  mu = xis[i];
89  x = 0.5 * ( x1 * ( 1 - mu ) + x2 * ( mu + 1 ) );
90  if( ( status = func( x, &sum, argList ) ) != nfu_Okay ) break;
91  *integral += sum * weights[i];
92  if( mu == 0 ) continue;
93  x = x1 + x2 - x;
94  if( ( status = func( x, &sum, argList ) ) != nfu_Okay ) break;
95  *integral += sum * weights[i];
96  }
97  }
98  *integral *= 0.5 * ( x2 - x1 );
99  return( status );
100 }
tuple x
Definition: test.py:50
static double sqrt_inv3
static constexpr double degree
Definition: G4SIunits.hh:144
enum nfu_status_e nfu_status
const G4int n
static struct nf_Legendre_GaussianQuadrature_degree GaussianQuadrature_degrees[nSets]

Here is the caller graph for this function:

double nf_Legendre_getCl ( nf_Legendre Legendre,
int  l,
nfu_status status 
)

Definition at line 119 of file nf_Legendre.cc.

119  {
120 
121  *status = nfu_Okay;
122  if( ( l < 0 ) || ( l > Legendre->maxOrder ) ) {
123  *status = nfu_badIndex;
124  return( 0. );
125  }
126  return( Legendre->Cls[l] );
127 }
double * Cls
Definition: nf_Legendre.h:26
int nf_Legendre_maxOrder ( nf_Legendre Legendre)

Definition at line 105 of file nf_Legendre.cc.

105  {
106 
107  return( Legendre->maxOrder );
108 }
nf_Legendre* nf_Legendre_new ( int  initialSize,
int  maxOrder,
double *  Cls,
nfu_status status 
)

Definition at line 23 of file nf_Legendre.cc.

23  {
24 
25  int l;
26  nf_Legendre *Legendre = (nf_Legendre *) nfu_malloc( sizeof( nf_Legendre ) );
27 
28  *status = nfu_mallocError;
29  if( Legendre == NULL ) return( NULL );
30  if( ( *status = nf_Legendre_setup( Legendre, initialSize, maxOrder ) ) != nfu_Okay ) {
31  nfu_free( Legendre );
32  return( NULL );
33  }
34  for( l = 0; l <= Legendre->maxOrder; l++ ) Legendre->Cls[l] = Cls[l];
35  return( Legendre );
36 }
nfu_status nf_Legendre_setup(nf_Legendre *nfL, int initialSize, int maxOrder)
Definition: nf_Legendre.cc:40
double * Cls
Definition: nf_Legendre.h:26
void * nfu_free(void *p)
void * nfu_malloc(size_t size)

Here is the call graph for this function:

Here is the caller graph for this function:

nfu_status nf_Legendre_normalize ( nf_Legendre Legendre)

Definition at line 146 of file nf_Legendre.cc.

146  {
147 
148  int l;
149  double norm;
150 
151  if( Legendre->maxOrder >= 0 ) {
152  if( ( norm = Legendre->Cls[0] ) == 0 ) return( nfu_divByZero );
153  for( l = 0; l <= Legendre->maxOrder; l++ ) Legendre->Cls[l] /= norm;
154  }
155  return( nfu_Okay );
156 }
double * Cls
Definition: nf_Legendre.h:26
double nf_Legendre_PofL_atMu ( int  l,
double  mu 
)

Definition at line 175 of file nf_Legendre.cc.

175  {
176 
177  int l_, twoL_plus1;
178  double Pl_minus1, Pl, Pl_plus1;
179 
180  if( l == 0 ) {
181  return( 1. ); }
182  else if( l == 1 ) {
183  return( mu ); }
184 /*
185  else if( l <= 9 ) {
186  double mu2 = mu * mu;
187  if ( l == 2 ) {
188  return( 1.5 * mu2 - 0.5 ); }
189  else if( l == 3 ) {
190  return( 2.5 * mu2 - 1.5 ) * mu; }
191  else if( l == 4 ) {
192  return( 4.375 * mu2 - 3.75 ) * mu2 + 0.375; }
193  else if( l == 5 ) {
194  return( ( 7.875 * mu2 - 8.75 ) * mu2 + 1.875 ) * mu; }
195  else if( l == 6 ) {
196  return( ( 14.4375 * mu2 - 19.6875 ) * mu2 + 6.5625 ) * mu2 - 0.3125; }
197  else if( l == 7 ) {
198  return( ( ( 26.8125 * mu2 - 43.3125 ) * mu2 + 19.6875 ) * mu2 - 2.1875 ) * mu; }
199  else if( l == 8 ) {
200  return( ( ( 50.2734375 * mu2 - 93.84375 ) * mu2 + 54.140625 ) * mu2 - 9.84375 ) * mu2 + 0.2734375; }
201  else {
202  return( ( ( ( 94.9609375 * mu2 - 201.09375 ) * mu2 + 140.765625 ) * mu2 - 36.09375 ) * mu2 + 2.4609375 ) * mu;
203  }
204  }
205 */
206 
207  Pl = 0.;
208  Pl_plus1 = 1.;
209  for( l_ = 0, twoL_plus1 = 1; l_ < l; l_++, twoL_plus1 += 2 ) {
210  Pl_minus1 = Pl;
211  Pl = Pl_plus1;
212  Pl_plus1 = ( twoL_plus1 * mu * Pl - l_ * Pl_minus1 ) / ( l_ + 1 );
213  }
214  return( Pl_plus1 );
215 }

Here is the caller graph for this function:

nfu_status nf_Legendre_reallocateCls ( nf_Legendre Legendre,
int  size,
int  forceSmallerResize 
)

Definition at line 77 of file nf_Legendre.cc.

77  {
78 
79  nfu_status status = nfu_Okay;
80 
82  if( size > ( nf_Legendre_maxMaxOrder + 1 ) ) size = nf_Legendre_maxMaxOrder + 1;
83  if( size != Legendre->allocated ) {
84  if( size > Legendre->allocated ) {
85  Legendre->Cls = (double *) nfu_realloc( size * sizeof( double ), Legendre->Cls ); }
86  else {
87  if( size < ( Legendre->maxOrder + 1 ) ) size = Legendre->maxOrder + 1;
88  if( ( Legendre->allocated > 2 * size ) || forceSmallerResize ) {
89  Legendre->Cls = (double *) nfu_realloc( size * sizeof( double ), Legendre->Cls ); }
90  else {
91  size = Legendre->allocated;
92  }
93  }
94  if( Legendre->Cls == NULL ) {
95  size = 0;
96  status = nfu_mallocError;
97  }
98  Legendre->allocated = size;
99  }
100  return( status );
101 }
double * Cls
Definition: nf_Legendre.h:26
#define nf_Legendre_maxMaxOrder
Definition: nf_Legendre.h:18
enum nfu_status_e nfu_status
void * nfu_realloc(size_t size, void *old)
#define nf_Legendre_minMaxOrder
Definition: nf_Legendre.h:17

Here is the call graph for this function:

Here is the caller graph for this function:

nfu_status nf_Legendre_release ( nf_Legendre nfL)

Definition at line 52 of file nf_Legendre.cc.

52  {
53 
54  if( Legendre->allocated > 0 ) nfu_free( Legendre->Cls );
55  memset( Legendre, 0, sizeof( nf_Legendre ) );
56  return( nfu_Okay );
57 }
void * nfu_free(void *p)

Here is the call graph for this function:

Here is the caller graph for this function:

nfu_status nf_Legendre_setCl ( nf_Legendre Legendre,
int  l,
double  Cl 
)

Definition at line 131 of file nf_Legendre.cc.

131  {
132 
133  nfu_status status;
134 
135  if( ( l < 0 ) || ( l > ( Legendre->maxOrder + 1 ) ) ) return( nfu_badIndex );
136  if( Legendre->allocated <= l ) {
137  if( ( status = nf_Legendre_reallocateCls( Legendre, l + nf_Legendre_sizeIncrement, 0 ) ) != nfu_Okay ) return( status );
138  }
139  if( l > Legendre->maxOrder ) Legendre->maxOrder = l;
140  Legendre->Cls[l] = Cl;
141  return( nfu_Okay );
142 }
#define nf_Legendre_sizeIncrement
Definition: nf_Legendre.h:19
double * Cls
Definition: nf_Legendre.h:26
enum nfu_status_e nfu_status
nfu_status nf_Legendre_reallocateCls(nf_Legendre *Legendre, int size, int forceSmallerResize)
Definition: nf_Legendre.cc:77

Here is the call graph for this function:

Here is the caller graph for this function:

nfu_status nf_Legendre_setup ( nf_Legendre nfL,
int  initialSize,
int  maxOrder 
)

Definition at line 40 of file nf_Legendre.cc.

40  {
41 
42  memset( Legendre, 0, sizeof( nf_Legendre ) );
43  if( maxOrder < 0 ) maxOrder = -1;
44  if( maxOrder > nf_Legendre_maxMaxOrder ) maxOrder = nf_Legendre_maxMaxOrder;
45  Legendre->maxOrder = maxOrder;
46  if( initialSize < ( maxOrder + 1 ) ) initialSize = maxOrder + 1;
47  return( nf_Legendre_reallocateCls( Legendre, initialSize, 0 ) );
48 }
#define nf_Legendre_maxMaxOrder
Definition: nf_Legendre.h:18
nfu_status nf_Legendre_reallocateCls(nf_Legendre *Legendre, int size, int forceSmallerResize)
Definition: nf_Legendre.cc:77

Here is the call graph for this function:

Here is the caller graph for this function:

ptwXYPoints* nf_Legendre_to_ptwXY ( nf_Legendre nfL,
double  accuracy,
int  biSectionMax,
int  checkForRoots,
nfu_status status 
)

Definition at line 219 of file nf_Legendre.cc.

219  {
220 
221  int i, n = 1;
222  double dx, xs[1000];
223  void *argList = (void *) Legendre;
224 
225  *status = nfu_Okay;
226  xs[0] = -1;
227  if( Legendre->maxOrder > 1 ) {
228  n = Legendre->maxOrder - 1;
229  if( n > 249 ) n = 249;
230  n = 4 * n + 1;
231  dx = 2. / n;
232  for( i = 1; i < n; i++ ) xs[i] = xs[i-1] + dx;
233  }
234  xs[n] = 1.;
235  return( ptwXY_createFromFunction( n + 1, xs, nf_Legendre_to_ptwXY2, (void *) argList, accuracy, checkForRoots, biSectionMax, status ) );
236 }
ptwXYPoints * ptwXY_createFromFunction(int n, double *xs, ptwXY_createFromFunction_callback func, void *argList, double accuracy, int checkForRoots, int biSectionMax, nfu_status *status)
Definition: ptwXY_misc.cc:40
const G4int n
static nfu_status nf_Legendre_to_ptwXY2(double mu, double *P, void *argList)
Definition: nf_Legendre.cc:240

Here is the call graph for this function: