OSHUN  beta
Arbitrary Order Spherical-Harmonic 1D-3P Vlasov-Fokker-Planck-Maxwell code
Data Structures | Functions
Algorithms Namespace Reference

Data Structures

class  AbstFunctor
 
class  Axis
 
class  AxisBundle
 
class  CAxis
 
class  LEAPs
 
class  LEAPv
 
class  PEFRL
 
class  RK2
 
class  RK3
 
class  RK4
 

Functions

template<typename T >
valarray< T > MakeAxis (const T min, const T max, const size_t N)
 
template<typename T >
valarray< T > MakeCAxis (const T min, const T max, const size_t N)
 
template<class T >
valarray< T > Legendre (const T x, const size_t Nl)
 
template<class T >
Array2D< T > Legendre (const T x, const size_t Nl, const size_t Nm)
 
template<class T >
T moment (const vector< T > q, const vector< T > x, const int p)
 
template<class T >
T moment (const vector< T > q, const valarray< T > x, const int p)
 
template<class T >
T moment (const valarray< T > q, const valarray< T > x, const int p)
 
template<class T >
T relativistic_invg_moment (const vector< T > q, const valarray< T > x, const int p)
 
template<class T >
T relativistic_invg_moment (const valarray< T > q, const valarray< T > x, const int p)
 
template<class T >
T relativistic_gamma_moment (const valarray< T > q, const valarray< T > x, const int p)
 

Function Documentation

◆ Legendre() [1/2]

template<class T >
valarray<T> Algorithms::Legendre ( const T  x,
const size_t  Nl 
)

Definition at line 227 of file lib-algorithms.h.

Referenced by Output_Data::PLegendre1D::PLegendre1D(), and Output_Data::PLegendre2D::PLegendre2D().

227  { // where Nl is to denote l0+1
228 
229  valarray<T> P_Legendre(0.0, Nl);
230  if ( abs(x) > 1 ) return P_Legendre;
231 
232  T r1, sqrtx = sqrt(1.0-x*x);
233 
234 // Initialization
235  P_Legendre[0] = 1.0;
236  P_Legendre[1] = x;
237 
238  for (size_t l(1); l < Nl - 1; ++l){
239  r1 = 1.0 / double(l + 1);
240  P_Legendre[l+1] = P_Legendre[l] * (x*(2.0*l+1.0) * r1) -
241  P_Legendre[l-1]*(double(l) * r1);
242  }
243 
244  return P_Legendre;
245  }
Here is the caller graph for this function:

◆ Legendre() [2/2]

template<class T >
Array2D<T> Algorithms::Legendre ( const T  x,
const size_t  Nl,
const size_t  Nm 
)

Definition at line 249 of file lib-algorithms.h.

249  {
250 
251 // Local variables
252  if (Nl < Nm) {
253  cout << "ERROR: " << "l0+1 = " << Nl << " < " << Nm << " = m0+1\n";
254  exit(1);
255  }
256 
257  Array2D<T> P_Legendre(Nl+1,Nm+1);
258 
259  if ( abs(x) > 1 ) {
260 // cout << "ERROR: " << "|" << x << "| > 1 is not a valid cosine\n";
261 // cout << "setting values to 0.0";
262  P_Legendre = 0.0;
263 
264  }
265  else
266  {
267 
268  T r1, sqrtx = sqrt(1.0-x*x), fact = 1.0;
269 
270  // Initialization
271  P_Legendre = 0.0;
272  P_Legendre(0,0) = 1.0;
273 
274  for (size_t l = 1; l < Nm+1; ++l){
275  // std::cout << "\n l=" << l;
276  P_Legendre(l,l) = - P_Legendre(l-1,l-1)*(fact*sqrtx);
277  fact += 2.0;
278  }
279  // std::cout << "\n ";
280 
281  for (size_t l = 0; l < ((Nm < Nl) ? Nm+1 : Nl+1); ++l){
282  // std::cout << "\n l=" << l;
283  P_Legendre(l+1,l) = P_Legendre(l,l)*(x*(2.0*l+1.0));
284  }
285  // std::cout << "\n ";
286 
287  for (size_t m = 0; m < Nm+1; ++m){
288  for (size_t l = m+1; l < Nl ; ++l){
289  // std::cout << "\n l=" << l << ", m=" << m;
290  r1 = 1.0 / double(l - m + 1);
291  P_Legendre(l+1,m) = P_Legendre(l,m) * (x*(2.0*l+1.0) * r1) -
292  P_Legendre(l-1,m)*(double(l+m) * r1);
293  }
294  }
295  }
296 
297  return P_Legendre;
298  }

◆ MakeAxis()

template<typename T >
valarray<T> Algorithms::MakeAxis ( const T  min,
const T  max,
const size_t  N 
)

Definition at line 25 of file lib-algorithms.h.

Referenced by Algorithms::Axis< T >::Axis(), Output_Data::p2p1x1_1D::p2p1x1_1D(), Output_Data::PLegendre1D::PLegendre1D(), and Output_Data::PLegendre2D::PLegendre2D().

25  {
26  valarray<T> v(N);
27  for (size_t i(0); i < N; ++i) {
28  v[i] = static_cast<T>(i);
29  }
30  v *= (max-min)/(static_cast<T>(N-1));
31  v += min;
32  return v;
33  }
Here is the caller graph for this function:

◆ MakeCAxis()

template<typename T >
valarray<T> Algorithms::MakeCAxis ( const T  min,
const T  max,
const size_t  N 
)

Definition at line 36 of file lib-algorithms.h.

Referenced by Export_Files::Header::axis(), Algorithms::CAxis< T >::CAxis(), DistFunc1D::getcurrent(), DistFunc1D::getdensity(), DistFunc1D::getpressure(), Output_Data::Output_Preprocessor_1D::Jx(), Output_Data::Output_Preprocessor_1D::Jy(), Output_Data::Output_Preprocessor_1D::Jz(), Output_Data::Output_Preprocessor_1D::n(), Output_Data::Output_Preprocessor_1D::Qx(), Output_Data::Output_Preprocessor_1D::Qy(), Output_Data::Output_Preprocessor_1D::Qz(), Output_Data::Output_Preprocessor_1D::T(), Output_Data::Output_Preprocessor_1D::vNx(), Output_Data::Output_Preprocessor_1D::vNy(), and Output_Data::Output_Preprocessor_1D::vNz().

36  {
37  valarray<T> v(N);
38  for (size_t i(0); i < N; ++i) {
39  v[i] = static_cast<T>(i);
40  }
41  v *= (max-min)/(static_cast<T>(N));
42  v += min+0.5*(max-min)/(static_cast<T>(N));
43  return v;
44  }
Here is the caller graph for this function:

◆ moment() [1/3]

template<class T >
T Algorithms::moment ( const vector< T >  q,
const vector< T >  x,
const int  p 
)

Definition at line 305 of file lib-algorithms.h.

Referenced by interspecies_f00_explicit_step::calculateintegrals(), DistFunc1D::getdensity(), DistFunc1D::getpressure(), Output_Data::Output_Preprocessor_1D::Jx(), Output_Data::Output_Preprocessor_1D::Jy(), Output_Data::Output_Preprocessor_1D::Jz(), Output_Data::Output_Preprocessor_1D::n(), Output_Data::Output_Preprocessor_1D::Qx(), Output_Data::Output_Preprocessor_1D::Qy(), Output_Data::Output_Preprocessor_1D::Qz(), Output_Data::Output_Preprocessor_1D::T(), Output_Data::Output_Preprocessor_1D::vNx(), Output_Data::Output_Preprocessor_1D::vNy(), and Output_Data::Output_Preprocessor_1D::vNz().

305  {
306 
307  T integral(0.0);
308 // TODO: Integral values up to the zeroth cell and above the last cell
309 // integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
310 // * (x[1]-x[0]);
311 // for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
312 // integral += q[i] * pow(x[i], p)
313 // * (x[i+1] - x[i-1]);
314 // }
315 // integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
316 // * (x[q.size()-1]-x[q.size()-2]);
317 // return integral*0.5;
318 
319  integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
320  * (x[1]-x[0]);
321  for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
322  integral += q[i] * pow(x[i], p)
323  * 0.5*(x[i+1] - x[i-1]);
324  }
325  integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
326  * (x[q.size()-1]-x[q.size()-2]);
327  return integral;
328  }
Here is the caller graph for this function:

◆ moment() [2/3]

template<class T >
T Algorithms::moment ( const vector< T >  q,
const valarray< T >  x,
const int  p 
)

Definition at line 331 of file lib-algorithms.h.

331  {
332 
333  T integral(0.0);
334 // TODO: Integral values up to the zeroth cell and above the last cell
335 // integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
336 // * (x[1]-x[0]);
337 // for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
338 // integral += q[i] * pow(x[i], p)
339 // * (x[i+1] - x[i-1]);
340 // }
341 // integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
342 // * (x[q.size()-1]-x[q.size()-2]);
343 // return integral*0.5;
344  integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
345  * (x[1]-x[0]);
346  for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
347  integral += q[i] * pow(x[i], p)
348  * 0.5*(x[i+1] - x[i-1]);
349  }
350  integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
351  * (x[q.size()-1]-x[q.size()-2]);
352  return integral;
353  }

◆ moment() [3/3]

template<class T >
T Algorithms::moment ( const valarray< T >  q,
const valarray< T >  x,
const int  p 
)

Definition at line 355 of file lib-algorithms.h.

355  {
356 
357  T integral(0.0);
358 // TODO: Integral values up to the zeroth cell and above the last cell
359 // integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
360 // * (x[1]-x[0]);
361 // for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
362 // integral += q[i] * pow(x[i], p)
363 // * (x[i+1] - x[i-1]);
364 // }
365 // integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
366 // * (x[q.size()-1]-x[q.size()-2]);
367 // return integral*0.5;
368  integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
369  * (x[1]-x[0]);
370  for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
371  integral += q[i] * pow(x[i], p)
372  * 0.5*(x[i+1] - x[i-1]);
373  }
374  integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
375  * (x[q.size()-1]-x[q.size()-2]);
376  return integral;
377  }

◆ relativistic_gamma_moment()

template<class T >
T Algorithms::relativistic_gamma_moment ( const valarray< T >  q,
const valarray< T >  x,
const int  p 
)

Definition at line 445 of file lib-algorithms.h.

445  {
446 
447  T integral(0.0);
448 // TODO: Integral values up to the zeroth cell and above the last cell
449 // integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
450 // * (x[1]-x[0])
451 // * sqrt(1.0+x[0]*x[0]);
452 // for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
453 // integral += q[i] * pow(x[i], p)
454 // * (x[i+1] - x[i-1])
455 // * sqrt(1.0+x[i]*x[i]);
456 // }
457 // integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
458 // * (x[q.size()-1]-x[q.size()-2])
459 // * sqrt(1.0+x[q.size-1]*x[q.size-1]);
460 // return integral*0.5;
461  integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
462  * (x[1]-x[0])
463  * sqrt(1.0+x[0]*x[0]);
464  for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
465  integral += q[i] * pow(x[i], p)
466  * 0.5*(x[i+1] - x[i-1])
467  * sqrt(1.0+x[i]*x[i]);
468  }
469  integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
470  * (x[q.size()-1]-x[q.size()-2])
471  * sqrt(1.0+x[q.size-1]*x[q.size-1]);
472  return integral;
473  }

◆ relativistic_invg_moment() [1/2]

template<class T >
T Algorithms::relativistic_invg_moment ( const vector< T >  q,
const valarray< T >  x,
const int  p 
)

Definition at line 382 of file lib-algorithms.h.

Referenced by DistFunc1D::getcurrent().

382  {
383 
384  T integral(0.0);
385 // TODO: Integral values up to the zeroth cell and above the last cell
386  // integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
387  // * (x[1]-x[0])
388  // / sqrt(1.0+x[0]*x[0]);
389  // for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
390  // integral += q[i] * pow(x[i], p)
391  // * (x[i+1] - x[i-1])
392  // / sqrt(1.0+x[i]*x[i]);
393  // }
394  // integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
395  // * (x[q.size()-1]-x[q.size()-2])
396  // / sqrt(1.0+x[q.size()-1]*x[q.size()-1]);
397  // return integral*0.5;
398 
399  integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
400  * (x[1]-x[0])
401  / sqrt(1.0+x[0]*x[0]);
402  for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
403  integral += q[i] * pow(x[i], p)
404  * 0.5*(x[i+1] - x[i-1])
405  / sqrt(1.0+x[i]*x[i]);
406  }
407  integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
408  * (x[q.size()-1]-x[q.size()-2])
409  / sqrt(1.0+x[q.size()-1]*x[q.size()-1]);
410  return integral;
411  }
Here is the caller graph for this function:

◆ relativistic_invg_moment() [2/2]

template<class T >
T Algorithms::relativistic_invg_moment ( const valarray< T >  q,
const valarray< T >  x,
const int  p 
)

Definition at line 414 of file lib-algorithms.h.

414  {
415 
416  T integral(0.0);
417 // TODO: Integral values up to the zeroth cell and above the last cell
418  // integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
419  // * (x[1]-x[0])
420  // / sqrt(1.0+x[0]*x[0]);
421  // for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
422  // integral += q[i] * pow(x[i], p)
423  // * (x[i+1] - x[i-1])
424  // / sqrt(1.0+x[i]*x[i]);
425  // }
426  // integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
427  // * (x[q.size()-1]-x[q.size()-2])
428  // / sqrt(1.0+x[q.size-1]*x[q.size-1]);
429  // return integral*0.5;
430  integral += q[0] * pow(x[0], p) // += Q_0*x_0^p * (x_1-x_0)
431  * (x[1]-x[0])
432  / sqrt(1.0+x[0]*x[0]);
433  for (size_t i(1); i < q.size()-1; ++i){ // += Q_i*x_i^p * (x_{i+1}-x_{i-1})
434  integral += q[i] * pow(x[i], p)
435  * 0.5*(x[i+1] - x[i-1])
436  / sqrt(1.0+x[i]*x[i]);
437  }
438  integral += q[q.size()-1] * pow(x[q.size()-1], p) // += Q_n*x_n^p * (x_{n}-x_{n-1})
439  * (x[q.size()-1]-x[q.size()-2])
440  / sqrt(1.0+x[q.size-1]*x[q.size-1]);
441  return integral;
442  }