The met_rocket Program with Emulated Real Numbers Instrumented for Run-time Trace

Calls to subroutines have been inserted automatically into the code to capture the results of all numerical computations. These are the routines trace_r4_data, trace_i4_data etc. shown below. On the first run, the results are captured to a binary file, the trace file. On subsequent runs, the trace file is read by the trace routines and the newly computed results are compared with the stored results. In the case of real and complex numbers, if the results differ by a criterion percentage the statement is reported as a possible breakdown in precision.

Note that the routines which capture the results have a second integer argument which is a unique identifier for the call-site. This is used to identify the statements where computation has broken down, and also to verify that the program flow has not changed.

!H!**************************************************************************** !H! File: ../fpt_output\met_rocket.f90 !H! Output by FPT 3.8-j Win32 On 30-SEP-15 At 15:07:56 Input files: !H! Main: j:\fpt\emulated_reals\met_rocket\fpt\met_rocket_precision_rtt.fsp !H! Current: j:\fpt\emulated_reals\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 !--------------------------------------------------^--------------------------- !!! FPT - 3769 Symbol or expression replaced by numeric value !------------------------------------------------------------------------------ , p_burn_time = em_real_k4( 0.210999999999999988E+01) & ! secs !------------------------------------------------------------^----------------- !!! FPT - 3769 Symbol or expression replaced by numeric value !------------------------------------------------------------------------------ , p_burn_thrust = em_real_k4( 0.129476200000000012E+06) & ! lb-ft !--------------------------------------------------------------^--------------- !!! FPT - 3769 Symbol or expression replaced by numeric value !------------------------------------------------------------------------------ , p_initial_mass = em_real_k4( 0.672000000000000028E+02) & ! lb !---------------------------------------------------------------^-------------- !!! FPT - 3769 Symbol or expression replaced by numeric value !------------------------------------------------------------------------------ , p_stage_2_mass = em_real_k4( 0.1E+02) & ! lb !----------------------------------------------^------------------------------- !!! FPT - 3769 Symbol or expression replaced by numeric value !------------------------------------------------------------------------------ , p_fuel_mass = em_real_k4( 0.376000000000000014E+02) & ! lb !------------------------------------------------------------^----------------- !!! FPT - 3769 Symbol or expression replaced by numeric value !------------------------------------------------------------------------------ , p_launch_vel = em_real_k4(190.0) & ! fps , frametime = em_real_k4( 0.100000000000000002E-02) & ! secs !----------------------------------------------------------^------------------- !!! FPT - 3769 Symbol or expression replaced by numeric value !------------------------------------------------------------------------------ ) ! 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 trace_start_sub_program(1) ! ! ***************************************************************************** ! OPEN(UNIT = 7,FILE = 'met_rocket.txt',STATUS = 'UNKNOWN') DO u = 6,7 CALL trace_i4_data(u,181) 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 CALL trace_r4_data(dtr,182) frame = 0 CALL trace_i4_data(frame,183) time = 0.0D0 CALL trace_r4_data(time,184) hdd = 0.0D0 CALL trace_r4_data(hdd,185) hd = 0.0D0 CALL trace_r4_data(hd,186) h = 0.0D0 CALL trace_r4_data(h,187) xdd = 0.0D0 CALL trace_r4_data(xdd,188) xd = 0.0D0 CALL trace_r4_data(xd,189) x = 0.0D0 CALL trace_r4_data(x,190) phdd = 0.0D0 CALL trace_r4_data(phdd,191) phd = 0.0D0 CALL trace_r4_data(phd,192) pxdd = 0.0D0 CALL trace_r4_data(pxdd,193) pxd = 0.0D0 CALL trace_r4_data(pxd,194) ! ! Un-comment lines beginning CCC to control drag parameters ! ========================================================= ! WRITE(6,'(''Launch angle (degrees): '',$)') READ(5,*)thetad CALL trace_r4_data(thetad,195) theta = thetad*dtr CALL trace_r4_data(theta,196) ! c_drag(1,1) = 0.0D0 CALL trace_r4_data(c_drag(1,1),197) c_drag(2,1) = 0.0008D0 CALL trace_r4_data(c_drag(2,1),198) c_drag(1,2) = 0.0D0 CALL trace_r4_data(c_drag(1,2),199) c_drag(2,2) = 0.00021D0 CALL trace_r4_data(c_drag(2,2),200) monitor_interval = 1000 CALL trace_i4_data(monitor_interval,201) ! !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 CALL trace_i4_data(frame,202) ! IF(time.LE.p_burn_time)THEN stage = 1 CALL trace_i4_data(stage,203) mass = p_initial_mass- & MAX(p_fuel_mass*(time/p_burn_time),0.0) CALL trace_r4_data(mass,204) thrust = p_burn_thrust CALL trace_r4_data(thrust,205) ELSE stage = 2 CALL trace_i4_data(stage,206) mass = p_stage_2_mass CALL trace_r4_data(mass,207) thrust = 0.0D0 CALL trace_r4_data(thrust,208) ENDIF ! ! Pressure in atmospheres atm_pressure = 10.0D0**(-h/50850.0D0) CALL trace_r4_data(atm_pressure,209) ! vel = SQRT(xd*xd+hd*hd) CALL trace_r4_data(vel,210) ! Assume drag is linear with pressure and quadratic with velocity drag = atm_pressure* & (c_drag(2,stage)*vel**2+c_drag(1,stage)*vel) CALL trace_r4_data(drag,211) ! xdd = ((thrust-drag)/mass)*COS(theta) CALL trace_r4_data(xdd,212) hdd = ((thrust-drag)/mass)*SIN(theta)-g CALL trace_r4_data(hdd,213) ! ! Trapezoidal integration - AB2 xd = xd+frametime*(2*xdd-pxdd) CALL trace_r4_data(xd,214) hd = hd+frametime*(2*hdd-phdd) CALL trace_r4_data(hd,215) pxdd = xdd CALL trace_r4_data(pxdd,216) phdd = hdd CALL trace_r4_data(phdd,217) x = x+frametime*(2*xd-pxd) CALL trace_r4_data(x,218) h = h+frametime*(2*hd-phd) CALL trace_r4_data(h,219) pxd = xd CALL trace_r4_data(pxd,220) phd = hd CALL trace_r4_data(phd,221) IF(vel.GT.p_launch_vel)THEN theta = ATAN2(hd,xd) CALL trace_r4_data(theta,222) ENDIF time = time+frametime CALL trace_r4_data(time,223) ! IF(frame.EQ.monitor_interval*(frame/monitor_interval))THEN DO u = 6,7 CALL trace_i4_data(u,224) 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 CALL trace_i4_data(u,225) 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) CALL trace_terminate(2) ! ! ***************************************************************************** ! END ! *****************************************************************