ENCMP 100 Computer Programming for Engineers © Joseph et al., 2020–21
Lab Assignment 5: Mercury’s Perihelion
Precession and General Relativity
According to Wikipedia, general relativity “is the geometric theory of gravitation published by Albert
Einstein in 1915 and is the current description of gravitation in modern physics.” In this lab assignment,
a student completes a MATLAB program to test with data an accurate prediction of Einstein’s theory,
namely the perihelion precession of Mercury. Mercury’s orbit around the Sun is not a stationary ellipse,
as Newton’s theory predicts when there are no other bodies. With Einstein’s theory, the relative angle
of Mercury’s perihelion (position nearest the Sun) varies by about 575.31 arcseconds per century.
Version 0: Get Started
Unzip the files in the V0GetStarted.zip file into your Current Folder. The horizons_results.txt
file, which you may open and review in a text editor, is a relatively large file of solar system data, called
ephemerides, downloaded from a NASA HORIZONS Web-Interface using the following parameters.
Parameter Setting
Ephemeris Type: VECTORS
Target Body: Mercury [199]
Coordinate Origin: Sun (body center) [500@10]
Time Span: Start=1800-01-01, Stop=2000-12-31, Step=1 d
Table Settings: quantities code=1; output units=KM-S; labels=NO;
CSV format=YES; object page=NO
Display/Output: download/save (plain text file)
Open the perihelion.m file in the Editor Window. The program consists of one script and seven local
functions in one M-file. While the file has a comment header, no function has a comment header. Run
the program. As it runs, it outputs a few lines of text to the Command Window to indicate progress. It
also creates a plot in a Figure Window and a file, horizons_results.tif, in the Current Folder.
Open the horizons_results.tif file outside of MATLAB, e.g., using Windows Photo Viewer. It is a
binary file of the plot shown in the Figure Window. Compare it to the horizons_results_v0.tif file
from the V0GetStarted.zip file. The two .tif files should look (and be otherwise) identical.
Apply your knowledge and skill to analyze the initial script and additional functions of the perihelion
M-file. Use the debugger to step through the program and to examine, in the base workspace, the data
cell array after each of the loaddata, locate, and select functions is invoked by the script.
Version 1: Load and Select
The implementations of str2cell and select simulate their functional requirements to help students
Get Started. Rewrite these functions to meet their implied functional requirements without changing
other code (or the data). First, replace any tilde (~) in a function header with a variable name.
ENCMP 100 Computer Programming for Engineers © Joseph et al., 2020–21
To properly load the data, the str2cell helper function must extract a numeric date (a scalar of class
'double'), e.g., 2378496.5, a short string date (a vector of class 'char'), e.g., '1800-Jan-01', and
3D numeric coordinates (a row vector of class 'double') from a long string, passed to it from the main
loaddata function. These three datums must be returned as cells of a single 1-by-3 cell array.
The select function must return an output cell array, data, with a subset of the rows in its input cell
array, data. If the second input argument is called ystep, selected rows have years whose remainder is
zero after conversion to a number and division by ystep (use mod). Months of the selected rows must
also match (at least) one string of the third input argument, a cell array of character vectors.
Unzip the one file, horizon_results_v1.tif, in the V1Load&Select.zip file into your Current
Folder. When you complete the Version 1 functional requirements, the output file (and associated plot)
of your perihelion program, horizons_results.tif, should look like this file when opened.
Write suitable comment headers, in your own words, for the following four functions: loaddata,
str2cell, locate, and select. A function's comment header should summarize its purpose, input
and output arguments, and side effects (such as Command Window input and output, file input and
output, and Figure Window plots). Complete the initial comment header of the program.
Submit your Version 1 solution via eClass by the Version 1 deadline. Submit the perihelion M-file
only. Before submission, test it in a folder where all other relevant files are as originally given.
Version 2: Refine and Save
With Version 1, not only is the perihelion precession inaccurate but also it is imprecise. The best fit slope
is far from 575.31 arcsec/cent and the precession computed from actual data looks almost random. To
address the imprecision, we refine data cells at the selected perihelia given more precise text files.
Unzip the V2Refine&Save.zip file into your Current Folder. Each horizons_results_YYYY-MMM-
DD.txt file, which you may open and review in a text editor, is a relatively small file downloaded from
NASA HORIZONS with the same parameters as before except for the Time Span parameter (see below).
This precise data is only given at half-century intervals (not at selected quarter-century intervals).
Additional Data (File Name) Web Interface (Time Span Setting)
mins_positions_1800-03-21.txt Start=1800-03-20, Stop=1800-03-22, Step=1 m
mins_positions_1850-01-28.txt Start=1850-01-27, Stop=1850-01-29, Step=1 m
mins_positions_1900-03-04.txt Start=1900-03-03, Stop=1900-03-05, Step=1 m
mins_positions_1950-01-11.txt Start=1950-01-10, Stop=1950-01-12, Step=1 m
mins_positions_2000-02-16.txt Start=2000-02-15, Stop=2000-02-17, Step=1 m
To refine the data, first add the statement “data = refine(data,'horizons_results');” to the
script immediately after the select function is invoked. Next, create a corresponding function stub in
the same perihelion M-file. Then, implement the function to meet its implied requirements without
modifying the data files or other code. When the implementation is sufficiently correct, the program
creates a horizons_results.tif file that looks like the horizons_results_v2.tif file.
ENCMP 100 Computer Programming for Engineers © Joseph et al., 2020–21
The refine function returns a more precise subset of its first input argument, data. The only rows of
the returned cell array, data, are the ones for which a file exists, whose name concatenates the second
input argument and text in a data cell. Reuse the existing loaddata function! If loaddata returns a
non-empty cell array, locate the row at which Mercury’s distance from the Sun is a minimum. Use the
min function with a second output argument. You may also use the cell2mat function.
Add the statement “savedata(data,'horizons_results')” to the script, immediately after the
makeplot statement. Add a savedata local function to output a comma-separated-value file, a text
file (supported by Microsoft Excel), that saves the first argument, data, in the format implied by
horizons_results_v2.csv. File extension aside, the name of the output file equals the second
argument. A sufficiently correct program outputs a horizons_results.csv file that equals the
horizons_results_v2.csv file. Use low-level file output (fopen, fprintf, and fclose)!
After refinement, the plot illustrates improved precision and accuracy. Moreover, the actual data used
by the plot is summarized in a file accessible outside MATLAB. Write suitable comment headers, in your
own words, for the following five functions: refine, makeplot, arcsec, annotate, and savedata.
Review the initial comment header of the entire program, especially reported percentages.
Submit your Version 2 solution via eClass by the Version 2 deadline. Submit the perihelion M-file
only. Before submission, test it in a folder where all other relevant files are as originally given.
Revision History
This document was created in 2021 by Dileepan Joseph and reviewed by Edward Tiong and Jason Myatt.
The concept and some code derive from a Dean’s Research Award project (2018) by Joseph and Ragur
Krishnan, an undergraduate student and winner (2017) of the ENCMP 100 Programming Contest.






























































































学霸联盟