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 

Complex issue with kind=1

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



Joined: 18 May 2012
Posts: 682
Location: Hamilton, Lanarkshire, Scotland.

PostPosted: Sun Apr 09, 2023 9:11 pm    Post subject: Complex issue with kind=1 Reply with quote

The following program returns the correct results when compiled with Win32.

When compiled with X64 the results are incorrect. The error occurs at line 8, in which there is the product of two complex numbers, the latter defined by two nested function calls.

Code:

module problem_mod
integer, parameter :: wp = kind(1.0)
contains
! Rotate phasor by specified angle (deg)
  complex(kind=wp) function crotphasordeg(phasor,angle)
  complex(kind=wp), intent(in) :: phasor
  real(kind=wp), intent(in) :: angle
    crotphasordeg = phasor*cang2phasor(sdeg2rad(angle))   ! ERROR OCCURS IN THIS EXPRESSION
!    crotphasordeg = cang2phasor(sdeg2rad(angle))          ! This alternative form gives the         
!    crotphasordeg = phasor*crotphasordeg                  ! correct result.
  end function crotphasordeg

! Define the unit phasor at the specified angle
  complex(kind=wp) function cang2phasor(angle)
  real(kind=wp), intent(in) :: angle
    cang2phasor = exp(cmplx(0.0,1.0,kind=wp)*angle)    ! e^(j*angle_rad)
  end function cang2phasor

! Convert supplied angle in deg to radians
  real(kind=wp) function sdeg2rad(angle)
  real(kind=wp), intent(in) :: angle
    sdeg2rad = angle*(4.0*atan(1.0))/180.0
  end function sdeg2rad
 
end module problem_mod

program p
use problem_mod
complex(kind=wp) z, z1
real(kind=wp) angle
z = cmplx(1.0,0.0,kind=wp)
angle = 30.0
z1 = crotphasordeg(z,angle)
print*, 'The phasor ', z, 'rotated by ', angle, 'degrees is', z1
angle = -30.0
z1 = crotphasordeg(z,angle)
print*, 'The phasor ', z, 'rotated by ', angle, 'degrees is', z1
print*
end program p


If the module parameter wp is changed to kind(1.d0) everything is as expected.

Evaluating the required values in a single line of code (kind=1), i.e. no function calls returns the correct results for both win32 and x64.
Code:
print*, '+30', cmplx(1.0,0.0)*exp(cmplx(0.0,1.0)*(+30.0*(4.0*atan(1.0))/180.0))
print*, '-30', cmplx(1.0,0.0)*exp(cmplx(0.0,1.0)*(-30.0*(4.0*atan(1.0))/180.0))
end
Back to top
View user's profile Send private message Visit poster's website
PaulLaidler
Site Admin


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

PostPosted: Mon Apr 10, 2023 7:13 am    Post subject: Reply with quote

Ken

Thank you for the bug report. It appears to have been already fixed and I am assuming that it is the same bug as that recently reported by mecej4 elsewhere on the Forum.

This fix will be in the next release of FTN95.
Back to top
View user's profile Send private message AIM Address
Kenneth_Smith



Joined: 18 May 2012
Posts: 682
Location: Hamilton, Lanarkshire, Scotland.

PostPosted: Mon Apr 10, 2023 10:38 am    Post subject: Reply with quote

Paul, thatís good news, i.e. you donít see the real part of the result going off to an very large number.

I didnít make the connection with mecej4ís post in the SLATEC thread, but looking at it this morning I can see some similarities. In particular, the error I was wrestling with yesterday is sensitive to the requested rotation angle.

Iíll park this one for now, pending the next FTN95 release.

Ken
Back to top
View user's profile Send private message Visit poster's website
Display posts from previous:   
Post new topic   Reply to topic    forums.silverfrost.com Forum Index -> Support 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