fpt Reference: EMULATE REAL ARITHMETIC

fpt and WinFPT Reference Manual - Command-line Commands

| SimCon Home | Ref Manual Home |

EMULATE REAL ARITHMETIC

Syntax:

[DO NOT] EMULATE REAL ARITHMETIC

Function:

This command was implemented to support the analyses described in Emulating Real Arithmetic in Fortran and in the publication Collins J. 2017, "Testing the Numerical Precisions Required to Execute Real World Programs", IJSEA Vol 8, No 2, March 2017.

fpt modifies the code so as to emulate all real and complex arithmetic. The changes made are:

A USE statement for the module module_emulate_real_arithmetic is inserted into every main compilation unit.
A call to the subroutine initialise_emulated_arithmetic is inserted as the first executable statement in every main program (Usually there is only one).
Every declaration of type REAL is changed to a declaration of TYPE (em_real_k4) or TYPE (em_real_k8). Similarly complex declarations are changed to declarations of TYPE (em_complex_k4) or TYPE (em_complex_k8) where the k4 or k8 indicates the kind in the original declaration. (Note that the current implementation supports only two real kinds, but additional kinds will be added if they are required in the future.)
Declarations of the form described above are inserted for all implicitly declared REAL and COMPLEX objects.
Real and complex literal numbers in PARAMETER statements, DATA statements and data initialisation in declaration statements are enclosed in type constructors, to convert them to the emulated data types.
Invocations of the intrinsic functions MAX and MIN with large numbers of arguments are sub-divided so that no more than four arguments occur in any one invocation. Literal values within the invocations are promoted to the real kind with the highest precision of the other arguments. This prevents a factorial explosion of the number of overloads required. The change affects very few statements.
Intrinsic functions which change the real kind of the result are re-programmed as two nested intrinsics (This change also affects very few statements).
The module module_emulate_real_arithmetic is added to the project. A version of this module is distributed in the fpt/examples or WinFPT\examples directory.

Note that no other changes are required to the code if, and only if, the emulated real and complex types have a single component of the same type and kind as the original real and complex objects. The first implementation of emulated arithmetic was used to modify the precision of the real and complex values. In this implementation, the emulated types satisfied this requirement. However, if the emulated types are required to have multiple components additional changes are made to the code in I/O statements, and to handle EQUIVALENCE relationships.

The Module module_emulate_real_arithmetic

This module contains:

Declarations of the four emulated types, em_real_k4, em_real_k8, em_complex_k4 and em_complex_k8.
Overloads for all assignment, arithmetic and relational operators for all combinations of these types with each other and with all possible real and integer kinds which may be written in the code. There are some hundreds of overloads.
Overloads for all intrinsic functions and subroutines which have real or complex arguments or results. Again there are many.
The routines required to carry out the desired experiments. In the first implementation, these routines modified the precisions of the input and output arguments for each operation so as to control the number of mantissa bits. Calls to the experiment routines are inserted in all of the emulation routines.

Two additional components were added to the code in the investigation of precision. They are module_em_signals.f90 and signal_underflow.f90 . These were used to trap underflow caused by reducing the number of mantissa bits. They do not form a necessary part of the emulation process.

Where to Use this Command

Operating system command line	Yes
Configuration file, config.fsp	No
Specification (fsp) files, *.fsp	Yes
Interactively, to FPT> prompt	No
Interactive command files	No
Embedded in the Fortran code	No

Default

Real arithmetic is not emulated by default.

Examples

The complete code of the module module_emulate_real_arithmetic.f90, set up to examine different numerical precisions, the code of example programs and the associated fpt scripts will be found in the fpt distribution in the directories:

examples/emulated_reals/
examples/met_rocket/
examples/ss/

The original code of the program met_rocket.f is shown below:

C ***************************************************************************** C met_rocket.for 13-Jul-06 Stephen Day Mandy et al C 18-Oct-06 John Collins C C ***************************************************************************** C 18-Oct-06 John Collins Removed the (very few) F90 features C ***************************************************************************** C PROGRAM met_rocket C IMPLICIT NONE C C ***************************************************************************** C REAL 1 g 1 ,p_burn_time 1 ,p_burn_thrust 1 ,p_initial_mass 1 ,p_stage_2_mass 1 ,p_fuel_mass 1 ,p_launch_vel 1 ,frametime PARAMETER ( 1 g = 32.2D0 ! fpsps 1 ,p_burn_time = 2.11D0 ! secs 1 ,p_burn_thrust = 4021.0D0*g ! lb-ft 1 ,p_initial_mass = 67.2D0 ! lb 1 ,p_stage_2_mass = 10.0D0 ! lb 1 ,p_fuel_mass = 37.6D0 ! lb 1 ,p_launch_vel = 190.0 ! fps 1 ,frametime = 0.001D0 ! secs 1 ) C INTEGER 1 frame 1 ,monitor_interval 1 ,stage ! Stage, before/after separation 1 ,u ! Logical unit number C REAL 1 dtr 1 ,time 1 ,xdd ! Down range acceleration, fpsps 1 ,xd ! Down range velocity, f/s 1 ,x ! Down range distance, ft 1 ,pxdd ! Last frame accelleration 1 ,pxd ! Last frame velocity 1 ,hdd ! Height acceleration, fpsps 1 ,hd ! Height velocity, fps 1 ,h ! Height, ft 1 ,phdd ! Last frame height acceleration 1 ,phd ! Last frame height velocity 1 ,theta ! Angle of flight to ground, radians 1 ,thetad ! Launch angle to ground, degrees 1 ,thrust 1 ,mass ! lb 1 ,drag 1 ,vel ! Direction of flight 1 ,atm_pressure ! Atmospheres 1 ,c_drag(2,2) ! Drag coeff, linear/quad, stage 1/2 C C ***************************************************************************** C OPEN (UNIT = 7,FILE = 'met_rocket.txt',STATUS = 'UNKNOWN') DO u = 6,7 WRITE (u,'(/,''Met Rocket Simulation'')') WRITE (u,'(''====================='',/)') ENDDO WRITE (6,'(''The program writes the output data to the file '' 1 ,''met_rocket.txt in the current'',/,''directory. '' 1 ,''Edit the program as marked to get interactice '' 1 ,''control of all of'',/,''the parameters'',/)') C C Initial conditions dtr = ATAN(1.0D0)/45.0D0 frame = 0 time = 0.0D0 hdd = 0.0D0 hd = 0.0D0 h = 0.0D0 xdd = 0.0D0 xd = 0.0D0 x = 0.0D0 phdd = 0.0D0 phd = 0.0D0 pxdd = 0.0D0 pxd = 0.0D0 C C Un-comment lines beginning CCC to control drag parameters C ========================================================= C WRITE (6,'(''Launch angle (degrees): '',$)') READ (5,*)thetad theta = thetad*dtr C c_drag(1,1) = 0.0D0 c_drag(2,1) = 0.0008D0 c_drag(1,2) = 0.0D0 c_drag(2,2) = 0.00021D0 monitor_interval = 1000 C CCC WRITE(6,'(''Drag stage 1 linear coeff: '',$)') CCC READ(5,*)c_drag(1,1) CCC WRITE(6,'(''Drag stage 1 quadratic coeff: '',$)') CCC READ(5,*)c_drag(2,1) CCC WRITE(6,'(''Drag stage 2 linear coeff: '',$)') CCC READ(5,*)c_drag(1,2) CCC WRITE(6,'(''Drag stage 2 quadratic coeff: '',$)') CCC READ(5,*)c_drag(2,2) CCC WRITE(6,'(''Monitor interval: '',$)') CCC READ(5,*)monitor_interval C WRITE (7,'(''Lanch angle:'',T32,F12.2)')thetad WRITE (7,'(''Drag coefficients'')') WRITE (7,'('' Stage 1 linear:'',T32,F12.6)')c_drag(1,1) WRITE (7,'('' Stage 1 quadratic:'',T32,F12.6)')c_drag(2,1) WRITE (7,'('' Stage 2 linear:'',T32,F12.6)')c_drag(1,2) WRITE (7,'('' Stage 2 quadratic:'',T32,F12.6,/)')c_drag(2,2) C C ***************************************************************************** C DO WHILE (theta .GE. 0.0D0) frame = frame+1 C IF (time .LE. p_burn_time) THEN stage = 1 mass = p_initial_mass- 1 MAX(p_fuel_mass*(time/p_burn_time),0.0) thrust = p_burn_thrust ELSE stage = 2 mass = p_stage_2_mass thrust = 0.0D0 ENDIF C C Pressure in atmospheres atm_pressure = 10.0D0**(-h/50850.0D0) C vel = SQRT(xd*xd+hd*hd) C Assume drag is linear with pressure and quadratic with velocity drag = atm_pressure* 1 (c_drag(2,stage)*vel**2+c_drag(1,stage)*vel) C xdd = ((thrust-drag)/mass)*COS(theta) hdd = ((thrust-drag)/mass)*SIN(theta)-g C C Trapezoidal integration - AB2 xd = xd+frametime*(2*xdd-pxdd) hd = hd+frametime*(2*hdd-phdd) pxdd = xdd phdd = hdd x = x+frametime*(2*xd-pxd) h = h+frametime*(2*hd-phd) pxd = xd phd = hd IF (vel .GT. p_launch_vel) THEN theta = ATAN2(hd,xd) ENDIF time = time+frametime C IF (frame .EQ. monitor_interval*(frame/monitor_interval)) THEN DO u = 6,7 WRITE (u,'(/,''Frame:'',I8,'' Time:'',F10.2 1 ,'' Theta:'',F10.2,'' Pressure:'',F10.4)') 1 frame,time,theta/dtr,atm_pressure WRITE (u,'('' h '' 1 ,'' hd '' 1 ,'' hdd '' 1 ,'' x '' 1 ,'' xd '' 1 ,'' xdd '' 1 )') WRITE (u,'(6F13.5)')h,hd,hdd,x,xd,xdd ENDDO CCC WRITE (6,'(''<CR> to continue'')') CCC READ (5,*) ENDIF C ENDDO C C Final values DO u = 6,7 WRITE (u,'(//,''Maximum altitude'')') WRITE (u,'(/,''Frame:'',I8,'' Time:'',F10.2 1 ,'' Theta:'',F10.2,'' Pressure:'',F10.4)') 1 frame,time,theta/dtr,atm_pressure WRITE (u,'('' h '' 1 ,'' hd '' 1 ,'' hdd '' 1 ,'' x '' 1 ,'' xd '' 1 ,'' xdd '' 1 )') WRITE (u,'(6F13.5)')h,hd,hdd,x,xd,xdd WRITE (u,'(/,''End of run'')') ENDDO CLOSE (7) C C ***************************************************************************** C END ! *****************************************************************

This code is modified automatically by the command EMULATE REAL ARITHMETIC to:

!H!**************************************************************************** !H! File: ../fpt_output\met_rocket.f90 !H! Output by FPT 4.0-c Win32 On 30-JUN-17 At 22:58:10 Input files: !H! Main: e:\winfpt\fpt\precision\met_rocket\fpt\met_rocket_precision.fsp !H! Current: e:\winfpt\fpt\precision\met_rocket\original_source\met_rocket.f !H! Licensee: SimCon: Development version. !H!**************************************************************************** ! ***************************************************************************** ! met_rocket.for 13-Jul-06 Stephen Day Mandy et al ! 18-Oct-06 John Collins ! ! ***************************************************************************** ! 18-Oct-06 John Collins Removed the (very few) F90 features ! ***************************************************************************** ! PROGRAM met_rocket ! ! USE module_emulate_real_arithmetic ! IMPLICIT NONE ! ! ***************************************************************************** ! TYPE (em_real_k4) & g & ,p_burn_time & ,p_burn_thrust & ,p_initial_mass & ,p_stage_2_mass & ,p_fuel_mass & ,p_launch_vel & ,frametime PARAMETER ( & g = em_real_k4( 0.322000000000000028E+02) &! fpsps ,p_burn_time = em_real_k4( 0.210999999999999988E+01) &! secs ,p_burn_thrust = em_real_k4( 0.129476200000000012E+06) &! lb-ft ,p_initial_mass = em_real_k4( 0.672000000000000028E+02) &! lb ,p_stage_2_mass = em_real_k4( 0.1E+02) &! lb ,p_fuel_mass = em_real_k4( 0.376000000000000014E+02) &! lb ,p_launch_vel = em_real_k4(190.0) &! fps ,frametime = em_real_k4( 0.100000000000000002E-02) &! secs ) ! INTEGER & frame & ,monitor_interval & ,stage & ! Stage, before/after separation ,u ! Logical unit number ! TYPE (em_real_k4) & dtr & ,time & ,xdd & ! Down range acceleration, fpsps ,xd & ! Down range velocity, f/s ,x & ! Down range distance, ft ,pxdd & ! Last frame accelleration ,pxd & ! Last frame velocity ,hdd & ! Height acceleration, fpsps ,hd & ! Height velocity, fps ,h & ! Height, ft ,phdd & ! Last frame height acceleration ,phd & ! Last frame height velocity ,theta & ! Angle of flight to ground, radians ,thetad & ! Launch angle to ground, degrees ,thrust & ,mass & ! lb ,drag & ,vel & ! Direction of flight ,atm_pressure & ! Atmospheres ,c_drag(2,2) ! Drag coeff, linear/quad, stage 1/2 ! ! ***************************************************************************** ! CALL initialise_emulated_arithmetic ! OPEN (UNIT = 7,FILE = 'met_rocket.txt',STATUS = 'UNKNOWN') DO u = 6,7 WRITE (u,'(/,''Met Rocket Simulation'')') WRITE (u,'(''====================='',/)') ENDDO WRITE (6,'(''The program writes the output data to the file '' & & ,''met_rocket.txt in the current'',/,''directory. '' & & ,''Edit the program as marked to get interactice '' & & ,''control of all of'',/,''the parameters'',/)') ! ! Initial conditions dtr = ATAN(1.0D0)/45.0D0 frame = 0 time = 0.0D0 hdd = 0.0D0 hd = 0.0D0 h = 0.0D0 xdd = 0.0D0 xd = 0.0D0 x = 0.0D0 phdd = 0.0D0 phd = 0.0D0 pxdd = 0.0D0 pxd = 0.0D0 ! ! Un-comment lines beginning CCC to control drag parameters ! ========================================================= ! WRITE (6,'(''Launch angle (degrees): '',$)') READ (5,*)thetad theta = thetad*dtr ! c_drag(1,1) = 0.0D0 c_drag(2,1) = 0.0008D0 c_drag(1,2) = 0.0D0 c_drag(2,2) = 0.00021D0 monitor_interval = 1000 ! !CC WRITE(6,'(''Drag stage 1 linear coeff: '',$)') !CC READ(5,*)c_drag(1,1) !CC WRITE(6,'(''Drag stage 1 quadratic coeff: '',$)') !CC READ(5,*)c_drag(2,1) !CC WRITE(6,'(''Drag stage 2 linear coeff: '',$)') !CC READ(5,*)c_drag(1,2) !CC WRITE(6,'(''Drag stage 2 quadratic coeff: '',$)') !CC READ(5,*)c_drag(2,2) !CC WRITE(6,'(''Monitor interval: '',$)') !CC READ(5,*)monitor_interval ! WRITE (7,'(''Lanch angle:'',T32,F12.2)')thetad WRITE (7,'(''Drag coefficients'')') WRITE (7,'('' Stage 1 linear:'',T32,F12.6)')c_drag(1,1) WRITE (7,'('' Stage 1 quadratic:'',T32,F12.6)')c_drag(2,1) WRITE (7,'('' Stage 2 linear:'',T32,F12.6)')c_drag(1,2) WRITE (7,'('' Stage 2 quadratic:'',T32,F12.6,/)')c_drag(2,2) ! ! ***************************************************************************** ! DO WHILE (theta .GE. 0.0D0) frame = frame+1 ! IF (time .LE. p_burn_time) THEN stage = 1 mass = p_initial_mass- & MAX(p_fuel_mass*(time/p_burn_time),0.0) thrust = p_burn_thrust ELSE stage = 2 mass = p_stage_2_mass thrust = 0.0D0 ENDIF ! ! Pressure in atmospheres atm_pressure = 10.0D0**(-h/50850.0D0) ! vel = SQRT(xd*xd+hd*hd) ! Assume drag is linear with pressure and quadratic with velocity drag = atm_pressure* & (c_drag(2,stage)*vel**2+c_drag(1,stage)*vel) ! xdd = ((thrust-drag)/mass)*COS(theta) hdd = ((thrust-drag)/mass)*SIN(theta)-g ! ! Trapezoidal integration - AB2 xd = xd+frametime*(2*xdd-pxdd) hd = hd+frametime*(2*hdd-phdd) pxdd = xdd phdd = hdd x = x+frametime*(2*xd-pxd) h = h+frametime*(2*hd-phd) pxd = xd phd = hd IF (vel .GT. p_launch_vel) THEN theta = ATAN2(hd,xd) ENDIF time = time+frametime ! IF (frame .EQ. monitor_interval*(frame/monitor_interval)) THEN DO u = 6,7 WRITE (u,'(/,''Frame:'',I8,'' Time:'',F10.2 & & ,'' Theta:'',F10.2,'' Pressure:'',F10.4)') & frame,time,theta/dtr,atm_pressure WRITE (u,'('' h '' & & ,'' hd '' & & ,'' hdd '' & & ,'' x '' & & ,'' xd '' & & ,'' xdd '' & & )') WRITE (u,'(6F13.5)')h,hd,hdd,x,xd,xdd ENDDO !CC WRITE (6,'(''<CR> to continue'')') !CC READ (5,*) ENDIF ! ENDDO ! ! Final values DO u = 6,7 WRITE (u,'(//,''Maximum altitude'')') WRITE (u,'(/,''Frame:'',I8,'' Time:'',F10.2 & & ,'' Theta:'',F10.2,'' Pressure:'',F10.4)') & frame,time,theta/dtr,atm_pressure WRITE (u,'('' h '' & & ,'' hd '' & & ,'' hdd '' & & ,'' x '' & & ,'' xd '' & & ,'' xdd '' & & )') WRITE (u,'(6F13.5)')h,hd,hdd,x,xd,xdd WRITE (u,'(/,''End of run'')') ENDDO CLOSE (7) ! ! ***************************************************************************** ! END ! *****************************************************************

The code has been reformatted to free format, but apart from this, note the very small number of changes which have been made.