forums.silverfrost.com Forum Index forums.silverfrost.com
Welcome to the Silverfrost forums
 
 FAQFAQ   SearchSearch   MemberlistMemberlist   UsergroupsUsergroups   RegisterRegister 
 ProfileProfile   Log in to check your private messagesLog in to check your private messages   Log inLog in 

Precision of intrinsic math functions

 
Post new topic   Reply to topic    forums.silverfrost.com Forum Index -> General
View previous topic :: View next topic  
Author Message
bplevenson



Joined: 02 Nov 2006
Posts: 21

PostPosted: Mon Nov 06, 2006 3:40 pm    Post subject: Precision of intrinsic math functions Reply with quote

I declare all my variables and constants as real(kind=2), which should make them double precision. But when I apply intrinsic math routines such as sin(x) or log(x), the answer is only coming out in single precision. I tried declaring sin, etc. as if they were real(kind=2) variables, and it compiled, but there was no noticeable result. What am I doing wrong? Is there some compiler switch or something I can turn on to get double-precision results? I want to use this compiler for scientific research and I need the precision.

-BPL
Back to top
View user's profile Send private message
JohnCampbell



Joined: 16 Feb 2006
Posts: 2554
Location: Sydney

PostPosted: Mon Nov 06, 2006 6:49 pm    Post subject: Precision of intrinsic math functions Reply with quote

Barton,

The following extract provides standard conforming definition of kind parameters for most precisions offererd in Salford.

Why not just declare "REAL*8" as I find this the most unambiguous representation of precision. It may be non-standard but it is better.

John

! Precision test
!
! Define integer & real precision options P min
integer*4, parameter :: i2 = selected_int_kind (4) ! 4 <5
integer*4, parameter :: i4 = selected_int_kind (Cool ! 9 5+
integer*4, parameter :: i8 = selected_int_kind (12) ! 18 10+
integer*4, parameter :: r4 = selected_real_kind (p=6) ! 6 <7
integer*4, parameter :: r8 = selected_real_kind (p=12) ! 15 7+
!
! NOTE THE CHANGE FOR 64 bit arithmetic for rq
!
integer*4, parameter :: rq = selected_real_kind (p=15) ! 18 or 33 16+
!
! Use high precision
REAL (kind=rq) :: k4, k8, kq, e ! for higher precision test
integer (kind=i4) :: i, p, next
!
print *, 'Integer Kind i2, i4, i8', i2, i4, i8
print *, 'Real Kind r4, r8, rq', r4, r8, rq
!
PRINT *, '1_i2 precision=', int(log10(real(huge(1_i2))))
PRINT *, '1_i4 precision=', int(log10(real(huge(1_i4))))
PRINT *, '1_i8 precision=', int(log10(real(huge(1_i8))))
PRINT *, '0.1_r4 precision=', precision (0.1_r4)
PRINT *, '0.1_r8 precision=', precision (0.1_r8)
PRINT *, '0.1_rq precision=', precision (0.1_rq)
!
Back to top
View user's profile Send private message
PaulLaidler
Site Admin


Joined: 21 Feb 2005
Posts: 7916
Location: Salford, UK

PostPosted: Tue Nov 07, 2006 1:16 am    Post subject: Precision of intrinsic math functions Reply with quote

Barton

I have to admit that I do not understand your problem.
Are you saying that the following code does not give the correct double precision answer?

REAL(2) x, y
x = 0.3_2
y = sin(x)

Back to top
View user's profile Send private message AIM Address
bplevenson



Joined: 02 Nov 2006
Posts: 21

PostPosted: Tue Nov 07, 2006 4:19 am    Post subject: Precision of intrinsic math functions Reply with quote

[small]John Campbell wrote:[/small]
Why not just declare "REAL*8" as I find this the most unambiguous representation of precision. It may be non-standard but it is better.

I don't have a problem with double-precision variables. My problem is that the results of applying intrinsic math functions to them comes back with single-precision values. It's sin, cos and tan that are returning single-precision values, not the variables.

-BPL
Back to top
View user's profile Send private message
bplevenson



Joined: 02 Nov 2006
Posts: 21

PostPosted: Tue Nov 07, 2006 5:01 am    Post subject: Precision of intrinsic math functions Reply with quote

[small]Paul Laidler wrote:[/small]
I have to admit that I do not understand your problem.
Are you saying that the following code does not give the correct double precision answer?

REAL(2) x, y
x = 0.3_2
y = sin(x)


Essentially, yes. I wrote a routine to input an angle in degrees, convert to radians, and take the sine, cosine and tangent. For pi I declared it as a parameter with 19 digits of precision, then divided that by 180.0 for the degrees-to-radians conversion factor. The answers I get for 45 degrees, 60 degrees and 360 degrees are all off starting in the 8th significant digit -- i.e., where the end of the first four bytes comes. I'll post my code below. Sorry for the over-elaborate text to number conversion logic.

[pre]!=====
! Program to show use of intrinsic and user-written functions.
!=====
Program Trig

real(kind=2), parameter :: pi = 3.141592653589793

integer :: dotplace ! Location of decimal point.
integer :: ivalue ! Value as an integer.
character(len=6) :: MyFormat ! Conversion format, string to number.
character(len=80) :: reply ! The reply itself.
integer :: rlen ! Length of reply in characters.
real(kind=2) :: value ! Mathematical value of reply.

!-----

! Explain how the program works:

write(*, *) 'This program shows the result of using trig'
write(*, *) 'functions on values you enter. Enter angles'
write(*, *) 'in degrees. Enter BYE or QUIT to stop.'
write(*, *) ! This prints a blank line.

! Start loop by getting user input:

do
write(*, 10)
10 format('Value or stop code: =>', $)
read(*, *) reply

rlen = len_trim(reply)
call upcase@(reply)

! If reply is a quit code, quit:

if (reply == 'BYE' .or. reply == 'QUIT') then
write (*, *) 'Thanks. Bye now.'
exit ! Quit program.

! If reply is a number, convert it appropriately:

else
dotplace = index(reply, '.')

if (dotplace == 0) then ! No decimal point?
write(MyFormat, 20) rlen ! Use integer format.
20 format('(I', I1, ')')

read(reply, MyFormat) ivalue ! Convert string to number.
value = real(ivalue)
else
write(MyFormat, 30) rlen ! Use real number format.
30 format('(G', I1, '.1)')

read(reply, MyFormat) value ! Convert string to number.
end if

! Finally, do trig operations on the resulting value:

Call TrigStuff(value) ! Call subroutine.
end if
end do

contains

!=====
! TrigStuff applies the trig functions to the user-entered value.
!=====
subroutine TrigStuff(x)
real(kind=2) :: x

real(kind=2), parameter :: pi = 3.141592653589793238
real(kind=2), parameter :: DegToRad = pi / 180.0

real(kind=2) :: radians ! Angle in radians.

!-----

radians = x * DegToRad ! Convert to radians.

write(*, 10) x, real(dsin(radians))
write(*, 20) x, dcos(radians)
write(*, 30) x, dtan(radians)

10 format('Sine ', F10.5, ' = ', F12.9)
20 format('Cosine ', F10.5, ' = ', F12.9)
30 format('Tangent ', F10.5, ' = ', F12.9)

end subroutine TrigStuff

End Program Trig[/pre]
Back to top
View user's profile Send private message
PaulLaidler
Site Admin


Joined: 21 Feb 2005
Posts: 7916
Location: Salford, UK

PostPosted: Tue Nov 07, 2006 5:43 am    Post subject: Precision of intrinsic math functions Reply with quote

Barton

I do not have time to check out your program in detail but I would use...

180.0_2 instead of 180.0
sin instead of dsin
cos instead of dcos
tan instead of dtan.

None of these should make any difference but the alternative is to check your code line by line, stepping with the debugger.
Back to top
View user's profile Send private message AIM Address
bplevenson



Joined: 02 Nov 2006
Posts: 21

PostPosted: Tue Nov 07, 2006 6:30 am    Post subject: Precision of intrinsic math functions Reply with quote

[small]Paul Laidler wrote:[/small]
I do not have time to check out your program in detail but I would use...

180.0_2 instead of 180.0
sin instead of dsin
cos instead of dcos
tan instead of dtan.

None of these should make any difference but the alternative is to check your code line by line, stepping with the debugger.


I tried it. No effect.

I fiddled with the output formats a little to bring out the problem more clearly. Here's some sample output:

-----------------------------------------------------------
This program shows the result of using trig
functions on values you enter. Enter angles
in degrees. Enter BYE or QUIT to stop.


Value or stop code: =>45.0

Sine 45.000 degrees = 0.707106796640858
Cosine 45.000 degrees = 0.707106765732237
Tangent 45.000 degrees = 1.000000043711391

Value or stop code: =>bye
Thanks. Bye now.

Press RETURN to close window . . .
-----------------------------------------------------------



-BPL
Back to top
View user's profile Send private message
bplevenson



Joined: 02 Nov 2006
Posts: 21

PostPosted: Tue Nov 07, 2006 6:34 am    Post subject: Precision of intrinsic math functions Reply with quote

Hold the phone!!!

Paul, you were right. It is a bug, but not with the intrinsics.

I declared pi as

real(kind=2), parameter :: 3.141592653589793238

when I made it

real(kind=2), parameter :: 3.141592653589793238_2

suddenly I started getting the right answers:


--------------------------------------------------------
This program shows the result of using trig
functions on values you enter. Enter angles
in degrees. Enter BYE or QUIT to stop.


Value or stop code: =>45

Sine 45.000 degrees = 0.707106781186547
Cosine 45.000 degrees = 0.707106781186548
Tangent 45.000 degrees = 1.000000000000000

Value or stop code: =>bye
Thanks. Bye now.

Press RETURN to close window . . .
--------------------------------------------------------


In other words, it didn't matter that I declared the constant as real(kind=2). It only became real(kind=2) when I added the _2 suffix to the numeric literal in the declaration.

They really shouldn't have real*4 as the default for constants when it could be real*8 or even real*10. Hope they fix this in the next version.

-Happy in Pittsburgh
Back to top
View user's profile Send private message
PaulLaidler
Site Admin


Joined: 21 Feb 2005
Posts: 7916
Location: Salford, UK

PostPosted: Tue Nov 07, 2006 11:41 am    Post subject: Precision of intrinsic math functions Reply with quote

Barton

Try using the compiler option DEFAULT_REAL_KIND 2. The old name was DREAL.
Back to top
View user's profile Send private message AIM Address
Display posts from previous:   
Post new topic   Reply to topic    forums.silverfrost.com Forum Index -> General All times are GMT + 1 Hour
Page 1 of 1

 
Jump to:  
You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot vote in polls in this forum


Powered by phpBB © 2001, 2005 phpBB Group