monofonIC/external/panphasia_ho

modules on COSMA7 

intel_comp/2018    fftw/3.3.9cosma7               
intel_mpi/2018     gsl/2.5            









 The code calls a function to generate the k-space modes for
a portion of the Panphasia field given an input descriptor.




 Should be called early before significant memory is allocated.  It
uses quite a bit of memory itself, but tidies up afterwards.

Has OpenMP -     -DUSE_OPENMP     in the makefile



 The routines support both single and double precision
calculations in two senses.  

The Fourier computations can be single or double precision

MACROs    FFTW_REAL/FFTW_COMPLEX  used to define 'Fourier' precision types
                                                    float or double.

The Panphasia coefficients can be single or double precision

MACROs   PAN_REAL/PAN_COMPLEX define the Panphasia precision - either
                                                    float or double.

 
To change the Fourier precision edit PAN_FFTW3.h - by default
single precision unless 'FOURIER_DOUBLE' is defined.

To change the Panphasia precision edit panphasia_functions.h and
single precision unless 'PAN_DOUBLE_PRECISION' is defined.








 Code description
-------------------

makefile

CODE
----


main.c    - demo program only

pan_mpi_routines.c   - contains MPI calls



high_order_panphasia_routines.c    -  serial - contains some OpenMP
uniform_rand_threefry4x64.c        -  serial - random generator and tests

Include files
--------------

panphasia_functions.h    

PAN_FFTW3.h           - MACROS for single/double precision FTs

pan_matrices_order6.h   - matrix coefficients for 6th order scheme.


threefry.h        -  Random generator 
array.h  + features array .h files 




Development notes:

----------------------------------------------------

14th April 2021



   Found a bug in the OpenMP version.  Different numbers
of threads led to a subset of Fourier modes having
different values. The precise differences changed
each time the code was run.

 Debugged by turning of OpenMP section by section.
The section which uses the spherical bessel functions
turned out to be responsible.

 The faulty version collapsed for 4 loops over
multipole,x,y,z.   Changing this to a loop
over multipoles, and collapsing 3 coordinate
loops solved the problem.

The variable index1 of the return field
does not depend on the multipole, while
index2 does.  Both index1 and index2 are
private.  This means the return array
(index 1) is updated several times.
Presumably as these updated occur
in parallel with the 4 loop collapsed
version the return array was being
corrupted sometimes.


15th April
------------

  This version supercedes version given to 
Oliver to add to MonofonIC clone.

 Main difference is additional OpenMP
statements and the ability to specify
in the descriptor that modes less
than of equal to some dimensionless
integer wavenumber squared are set
to the mean power.

 Tested output on 1 core - with/without
OpenMP.  Not tested with more than
1 MPI rank.