Hi all.
The Cactus documentation shows a set of C routines CCTK_CmplxXXXXXX that presumably preserve precision between CCTK_REAL and CCTK_COMPLEX types. Is there a "right" way to handle this in Cactus with Fortran (that isn't "recode in C")?
Context: I have some Fortran code in a Cactus thorn that works with real & complex fields. I've noticed that when the numbers are very large in magnitude, taking a simple real(z) drops precision --- a finite but large real part ends up as +/- Inf instead. I suspect it's dropping to single precision at some stage. Switching to "dreal" doesn't help, and avoiding the use of CCTK_COMPLEX altogether --- while possible --- leads to longer, less readable code.
Thanks,
Bernard
Bernard
In Fortran, real(z) (where z is complex) does maintain precision. The problem is likely caused by your using "-Ofast", "-ffast-math", or an Intel compiler that uses the equivalent by default (so you have to explicitly switch it off). With such an optimization, e.g. abs(z) is implemented as sqrt(real(z)**2 + aimag(z)**2), which overflows the intermediate squares. Without -ffast-math, compilers use a more expensive implementation that does not overflow.
In your case, I'd assume that the problem is no the call to real, but other complex or real functions. Can you show us a larger part of the code, and the compiler vendor, version, and options you are using?
The C functions CCTK_Cmplx* are deprecated. These days, C compilers implement complex numbers as defined by the standard, and you'd use the respective C standard function instead.
-erik
On Sun, Jul 19, 2015 at 2:16 PM, Bernard Kelly physicsbeany@gmail.com wrote:
Hi all.
The Cactus documentation shows a set of C routines CCTK_CmplxXXXXXX that presumably preserve precision between CCTK_REAL and CCTK_COMPLEX types. Is there a "right" way to handle this in Cactus with Fortran (that isn't "recode in C")?
Context: I have some Fortran code in a Cactus thorn that works with real & complex fields. I've noticed that when the numbers are very large in magnitude, taking a simple real(z) drops precision --- a finite but large real part ends up as +/- Inf instead. I suspect it's dropping to single precision at some stage. Switching to "dreal" doesn't help, and avoiding the use of CCTK_COMPLEX altogether --- while possible --- leads to longer, less readable code.
Thanks,
Bernard _______________________________________________ Users mailing list Users@einsteintoolkit.org http://lists.einsteintoolkit.org/mailman/listinfo/users
Thanks, Erik.
I wasn't taking any magnitudes, but it looks like my problem was just blindly using the "cmplx" function to form the complex variable without a KIND argument, and it defaults to KIND=1 (i.e. single-precision) even if the real & imaginary components enter as double-precision.
For now, I've switched to the GNU "complex" method, which seems to have fixed my problem on my Macbook.
Bernard
On 19 July 2015 at 14:52, Erik Schnetter schnetter@cct.lsu.edu wrote:
Bernard
In Fortran, real(z) (where z is complex) does maintain precision. The problem is likely caused by your using "-Ofast", "-ffast-math", or an Intel compiler that uses the equivalent by default (so you have to explicitly switch it off). With such an optimization, e.g. abs(z) is implemented as sqrt(real(z)**2 + aimag(z)**2), which overflows the intermediate squares. Without -ffast-math, compilers use a more expensive implementation that does not overflow.
In your case, I'd assume that the problem is no the call to real, but other complex or real functions. Can you show us a larger part of the code, and the compiler vendor, version, and options you are using?
The C functions CCTK_Cmplx* are deprecated. These days, C compilers implement complex numbers as defined by the standard, and you'd use the respective C standard function instead.
-erik
On Sun, Jul 19, 2015 at 2:16 PM, Bernard Kelly physicsbeany@gmail.com wrote:
Hi all.
The Cactus documentation shows a set of C routines CCTK_CmplxXXXXXX that presumably preserve precision between CCTK_REAL and CCTK_COMPLEX types. Is there a "right" way to handle this in Cactus with Fortran (that isn't "recode in C")?
Context: I have some Fortran code in a Cactus thorn that works with real & complex fields. I've noticed that when the numbers are very large in magnitude, taking a simple real(z) drops precision --- a finite but large real part ends up as +/- Inf instead. I suspect it's dropping to single precision at some stage. Switching to "dreal" doesn't help, and avoiding the use of CCTK_COMPLEX altogether --- while possible --- leads to longer, less readable code.
Thanks,
Bernard _______________________________________________ Users mailing list Users@einsteintoolkit.org http://lists.einsteintoolkit.org/mailman/listinfo/users
-- Erik Schnetter schnetter@cct.lsu.edu http://www.perimeterinstitute.ca/personal/eschnetter/
Bernard
I usually define some constants in each routine or module:
CCTK_REAL, parameter :: zero = 0 integer, parameter :: rkind = kind(zero)
and then I can write "cmplx(x, kind=rk)".
-erik
On Mon, Jul 20, 2015 at 12:09 PM, Bernard Kelly physicsbeany@gmail.com wrote:
Thanks, Erik.
I wasn't taking any magnitudes, but it looks like my problem was just blindly using the "cmplx" function to form the complex variable without a KIND argument, and it defaults to KIND=1 (i.e. single-precision) even if the real & imaginary components enter as double-precision.
For now, I've switched to the GNU "complex" method, which seems to have fixed my problem on my Macbook.
Bernard
On 19 July 2015 at 14:52, Erik Schnetter schnetter@cct.lsu.edu wrote:
Bernard
In Fortran, real(z) (where z is complex) does maintain precision. The problem is likely caused by your using "-Ofast", "-ffast-math", or an
Intel
compiler that uses the equivalent by default (so you have to explicitly switch it off). With such an optimization, e.g. abs(z) is implemented as sqrt(real(z)**2 + aimag(z)**2), which overflows the intermediate squares. Without -ffast-math, compilers use a more expensive implementation that
does
not overflow.
In your case, I'd assume that the problem is no the call to real, but
other
complex or real functions. Can you show us a larger part of the code, and the compiler vendor, version, and options you are using?
The C functions CCTK_Cmplx* are deprecated. These days, C compilers implement complex numbers as defined by the standard, and you'd use the respective C standard function instead.
-erik
On Sun, Jul 19, 2015 at 2:16 PM, Bernard Kelly physicsbeany@gmail.com wrote:
Hi all.
The Cactus documentation shows a set of C routines CCTK_CmplxXXXXXX that presumably preserve precision between CCTK_REAL and CCTK_COMPLEX types. Is there a "right" way to handle this in Cactus with Fortran (that isn't "recode in C")?
Context: I have some Fortran code in a Cactus thorn that works with real & complex fields. I've noticed that when the numbers are very large in magnitude, taking a simple real(z) drops precision --- a finite but large real part ends up as +/- Inf instead. I suspect it's dropping to single precision at some stage. Switching to "dreal" doesn't help, and avoiding the use of CCTK_COMPLEX altogether --- while possible --- leads to longer, less readable code.
Thanks,
Bernard _______________________________________________ Users mailing list Users@einsteintoolkit.org http://lists.einsteintoolkit.org/mailman/listinfo/users
-- Erik Schnetter schnetter@cct.lsu.edu http://www.perimeterinstitute.ca/personal/eschnetter/
users@lists.einsteintoolkit.org