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 <[email protected]> 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 > [email protected] > http://lists.einsteintoolkit.org/mailman/listinfo/users > -- Erik Schnetter <[email protected]> http://www.perimeterinstitute.ca/personal/eschnetter/
_______________________________________________ Users mailing list [email protected] http://lists.einsteintoolkit.org/mailman/listinfo/users
