
forums.silverfrost.com Welcome to the Silverfrost forums

View previous topic :: View next topic 
Author 
Message 
DanRRight
Joined: 10 Mar 2008 Posts: 2867 Location: South Pole, Antarctica

Posted: Sun Sep 15, 2024 4:11 am Post subject: Any MATLAB Pros are here? 


Can you help to convert this Fortran code to Matlab?
Code:  integer,parameter::k=selected_int_kind(18)
integer(k) v
integer U,L
111 READ(*,*) v
U = ishft(v, 24)
L = iand(v, Z'000000ffffff')
print*,L,U
! getting v back
v = int(L,k) + ishft(int(U,k),24)
print*,v
goto 111
end 
We discussed this code a year ago, see https://forums.silverfrost.com/viewtopic.php?t=4840.
This code reads from screen 64bit number v you put from the keyboard and splits it into two 32bit single precision numbers U and L. Then it combines two numbers into v back to check if all went OK.
This code takes only first 24bits of the v number to restrict U and L to values below approximately 20 million, because if you convert more and then assign these U and L integer numbers to real single precision 32bit numbers (and we have to do that), all further bits will be lost (mantissa of real number is just 24 bits or so). This code with this tricky way can split v numbers up to 10^15 instead of just 20 million or so 

Back to top 


DanRRight
Joined: 10 Mar 2008 Posts: 2867 Location: South Pole, Antarctica

Posted: Mon Sep 16, 2024 3:45 pm Post subject: 


Seems there was not a Matlab but the FTN95 problem. We wrote Matlab code but it produced different result for U. I just tried Gfortran and also one additional Fortran compiler and both also "agree to disagree" with FTN95. If we try any v number less than 16 million FTN95 gives U=1 but others having U=0. All generate the same correct value L. Fun is that FTN95 somehow still recovers v value back correctly, which means it has double problem ? Or it is all other compilers which have a problem? 

Back to top 


PaulLaidler Site Admin
Joined: 21 Feb 2005 Posts: 8019 Location: Salford, UK

Posted: Mon Sep 16, 2024 4:16 pm Post subject: 


I don't have time to look at this now but I seem to recall that the last time we looked at this was http://forums.silverfrost.com/viewtopic.php?t=4840&highlight=ishft.
The Z constants required a kind specifier. There may be differences in the required Z form between one compiler and another. 

Back to top 


DanRRight
Joined: 10 Mar 2008 Posts: 2867 Location: South Pole, Antarctica

Posted: Tue Sep 17, 2024 6:22 am Post subject: 


That was good advice, thanks.
By the way Matlab simply refused to compile in this case giving an error. And Gfortran did not complain but gave correct result without me specifying 64bit Z number not as simply Znumber as Znumber_k.
May be it is still worth for FTN95 to make in this specific place an error because this trick is so obscure that only 1 per million people will find the reason of discrepancy. Problem also was that combining back the splitted numbers to check if splitting ended fine, was giving correct value  double error with correct income. Pure devilry 

Back to top 


JohnCampbell
Joined: 16 Feb 2006 Posts: 2587 Location: Sydney

Posted: Wed Sep 18, 2024 6:05 am Post subject: 


Dan,
I think the only problem is with iand, so I included a corrected I*8 mask.
Let me know if I have misunderstood.
I modified your code, to enforce a 64bit integer mask, "mask24".
I tested for a number of values "v", but all worked with your 32bit mask.
You didn't provide values for v. Do you have any values where v fails ?
Actually, any values larger than 2^55 will fail (U is ve), but any others ?
Why use 24 bit split ?
This is my expanded code to identify any problems is:
Code:  ! This works for all values smaller that 2^(24+31),
! ie no overflow on I84 U
integer, parameter :: i8=selected_int_kind(18)
integer(i8), parameter :: two = 2
integer(i8), parameter :: mask24 = two**241
integer(i8) :: vlist(6) = [ 16777216, 16777215, 134217728, 134217727, 1234567890, 1 ]
integer(i8) :: v
integer :: U,L, k
write (*,11) mask24, mask24, mask24
11 format ( 'Mask ',z20,i20,b33)
! vlist(5) = two ** 34  3 ! larger value
vlist(5) = two ** 58  3 ! larger value that will fail as longer than 55 bits
!111 READ(*,*) v
do k = 1,size(vlist)
! Old approach
v = vlist(k)
U = ishft ( v, 24 )
L = iand ( v, Z'000000ffffff' )
! print*,L,U
! getting v back
v = int(L,i8) + ishft(int(U,i8),24)
! print*,v
print*,'old',U,L,v,vlist(k)
! New approach with I*8 : mask24
v = vlist(k)
U = ishft ( v, 24 ) ! problem if U ve
L = iand ( v, mask24 )
! print*,L,U
! getting v back
v = int(L,i8) + ishft(int(U,i8),24)
! print*,v
if ( v == vlist(k) ) then
print*,'new',U,L,v,vlist(k)
else
print*,'new',U,L,v,vlist(k),' ERROR'
end if
end do
end 


Back to top 


DanRRight
Joined: 10 Mar 2008 Posts: 2867 Location: South Pole, Antarctica

Posted: Sun Sep 29, 2024 7:37 pm Post subject: 


24 bit split used because after V splitted to two integer numbers U and L they will be assigned to real*4 numbers RU and RL and saved to media (that was done for speed considerations, RU and RL will be a part of some huge multiGB or TB size arrays) . When file will be loaded again these RU and RL will be converted to V preserving all its first digits up to something like 10^15 value. This is not full Integer*8 but I hope that the 640...sorry I mean 10^15 64bit numbers will fit to everyone. How many Petabytes is it? 

Back to top 




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
