Degenerate Conic

Algorithms • Modern Fortran Programming • Orbital Mechanics

Sep 13, 2018

Intel Fortran Compiler 19.0

intel-logo-vector

Intel has just released v19 of the Intel Fortran Compiler (part of Intel Parallel Studio XE 2019). The new release adds some new features from the upcoming Fortran 2018 standard:

  • Coarray events
  • Intrinsic function coshape
  • Default accessibility for entities accessed from a module
  • Import Enhancements
  • All standard procedures in ISO_C_BINDING other than C_F_POINTER are now PURE

Presumably, it also include some bug fixes. Intel used to provide a list of bug fixes, but they seem to have stopping doing that for some reason.

See also

Aug 05, 2018

Circular Restricted Three-Body Problem

halo-300x235

Libration Points in the Earth-Moon System. From Farquhar, 1971.

The Circular Restricted Three-Body Problem (CR3BP) is a continuing source of delight for astrodynamicists. If you think that, after 200 years, everything that could ever be written about this subject has already been written, you would be wrong. The basic premise is to describe the motion of a spacecraft in the vicinity of two celestial bodies that are in circular orbits about their common barycenter. The "restricted" part means that the spacecraft does not affect the motion of the other two bodies (i.e., its mass is negligible).

The normalized equations of motion for the CR3BP are very simple (see the Fortran Astrodynamics Toolkit):

subroutine crtbp_derivs(mu,x,dx)

implicit none

real(wp),intent(in) :: mu
!! CRTBP parameter
real(wp),dimension(6),intent(in) :: x
!! normalized state [r,v]
real(wp),dimension(6),intent(out) :: dx
!! normalized state derivative [rdot,vdot]

real(wp),dimension(3) :: r1,r2,rb1,rb2,r,v,g
real(wp) :: r13,r23,omm,c1,c2

r = x(1:3) ! position
v = x(4:6) ! velocity
omm = one - mu
rb1 = [-mu,zero,zero] ! location of body 1
rb2 = [omm,zero,zero] ! location of body 2
r1 = r - rb1 ! body1 -> sc vector
r2 = r - rb2 ! body2 -> sc vector
r13 = norm2(r1)**3
r23 = norm2(r2)**3
c1 = omm/r13
c2 = mu/r23

!normalized gravity from both bodies:
g(1) = -c1*(r(1) + mu) - c2*(r(1)-one+mu)
g(2) = -c1*r(2) - c2*r(2)
g(3) = -c1*r(3) - c2*r(3)

! derivative of x:
dx(1:3) = v ! rdot
dx(4) = two*v(2) + r(1) + g(1) ! vdot
dx(5) = -two*v(1) + r(2) + g(2) !
dx(6) = g(3) !

end subroutine crtbp_derivs

dros2

There are many types of interesting periodic orbits that exist in this system. I've mentioned Distant Retrograde Orbits (DROs) here before. Probably the most well-known CR3BP periodic orbits are called "halo orbits". Halo orbits can be computed in a similar manner as DROs (see previous post). For halos, we can get a good initial guess using an analytic approximation, and then use this guess to target the real thing. See the halo orbit module in the Fortran Astrodynamics Toolkit for an example of this approach. We can also use a continuation type method to generate families of orbits (say, for the L2 libration point, we have a range of halo orbits from near-planars that "orbit" L2, to near-rectilinears with close approaches to the secondary body). In a realistic force model (with the real ephemeris of the Earth, Moon, and Sun), the orbits generated using CR3BP assumptions are not periodic, but we can use the CR3BP solution as an initial guess in a multiple-shooting type approach with a numerical nonlinear equation solver to generate quasi-periodic multi-rev solutions. See the references below for details.

I generated the following image of a bunch of planar periodic CR3BP orbits using a program incorporating the Fortran Astrodynamics Toolkit and OpenFrames (an interface to the open source OpenSceneGraph library that can be used to visualize spacecraft trajectories). The initial states for the various orbits were taken from this database. Eventually, I will make this program open source when it is finished. There's also an online CR3BP propagator here.

screenshot5

See also

Jul 14, 2018

CMLIB

nbs

The National Bureau of Standards (NBS) Core Math Library (CMLIB) is a "collection of high-quality, easily transportable Fortran subroutine sublibraries solving standard problems in many areas of mathematics and statistics". It was written in FORTRAN 77 and is available from a NIST FTP link. The first version (1.0) was released in March 1986, and the last update (3.0) was in 1988. The software is no longer maintained, and I assume one day the link will probably disappear. It was compiled mostly from other externally-available libraries available at the time (including BLAS, EISPACK, FISHPAK, FNLIB, FFTPACK, LINPACK, and QUADPACK) but there are also some original codes. My Bspline-Fortran package is a modernized update and extension of the DTENSBS routines from CMLIB. The core spline interpolation routines (slightly updated) still work great after more than 30 years.

The National Bureau of Standards (NBS) was renamed the National Institute of Standards and Technology (NIST) in August 1988 (a few months after CMLIB was last updated). The nice thing about these old government-produced codes is that they are public domain in the United States. Indeed any software written by a U.S. government employee (but not a government contractor) as part of their official duties is automatically public domain. I generally prefer permissive software licenses, and you can't get more permissive than that.

See also

Jun 06, 2018

Fortran + JSON + Python

python_json_fortran

There are various ways to interface Fortran and Python. However, many of the online examples on this subject seem to be content with a 40 year old variant of the language, and do not tend to focus on modern Fortran. One particular thing that is hard to figure out is how to send arbitrary or dynamic data from Python to Fortran (and vice versa)? This could be variable-sized arrays, variable-sized strings (or even a variable number of variable sized strings), or just any random data structure that contains different types of variables that aren't necessarily known at compile time. It is not easy to find out how to do this, but one option presented here is to use JSON as an intermediate format.

Sending a dict from Python to Fortran

Say we want to communicate an arbitrary Python dictionary structure to Fortran. We can do this fairly easily by dumping the dict to a JSON string, sending the string to Fortran, and then parsing the string in Fortran using JSON-Fortran. Then we have access to all the data in the dictionary in the Fortran code. An example of this is presented here.

First we need a function python_str_to_fortran() to convert a Python dict object into a c_char_p string that can be sent via the C interface to a Fortran routine:

from ctypes import *
import json

def python_str_to_fortran(s):
    return cast(create_string_buffer(s.encode()),c_char_p)

def python_dict_to_fortran(d):
    return python_str_to_fortran(json.dumps(d))

Note that we are using Python's ctypes module here to create a c_char_p string that can be passed to Fortran. I admit, this bit of code is the result of some trial and error to figure out what would work (the documentation is less than helpful and of course doesn't mention Fortran at all). Ctypes requires that the Fortran code be compiled as a shared library (a .DLL on Windows, or .so on Linux). We can then define the interfaces to the routines in the DLL we want to call from Python. As long as we get the types of the arguments right (which can be tricky for strings as will be shown) it works great. Consider this example:

# create a dict:
a = {'Generated in Python': True,
     'scalar': 1,
     'vector': [1,2,3],
     'string': 'hello'}

cp = python_dict_to_fortran(a) # creates the c_char_p

dll = CDLL ('lib/json.so') # load the shared library

# define the interface to one of the Fortran routines:
send_json_to_fortran = dll.send_json_to_fortran
send_json_to_fortran.argtypes = [POINTER(c_char_p)]
send_json_to_fortran.restype = None

send_json_to_fortran(byref(cp)) # send the string to fortran

Here we are creating a dictionary, converting it into a string, loading our Fortran DLL, defining the interface to the Fortran function we want to call (send_json_to_fortran()), and then calling it. On the Fortran side, we first define a routine to convert the input string into a JSON structure (using JSON-Fortran):

subroutine c_ptr_to_json(cp,json)

type(c_ptr),intent(in) :: cp
!! a `c_char_p` from python containing a JSON string.
type(json_file),intent(inout) :: json
!! the JSON data structure

character(len=:),allocatable :: fstr
!! string containing the JSON data

if (c_associated(cp)) then
    fstr = c_ptr_to_f_string(cp)
    call json%load_from_string(fstr) ! parse JSON
    deallocate(fstr)
end if

end subroutine c_ptr_to_json

Note that we are using the following routine to convert the Python c_char_p string into a normal Fortran string:

function c_ptr_to_f_string(cp) result(fstr)

type(c_ptr),intent(in) :: cp !! a `c_char_p` from python
character(len=:),allocatable :: fstr
!! the corresponding fortran string

integer :: ilen !! string length

ilen = strlen(cp) ! C library function

block
    ! convert the C string to a Fortran string
    character(kind=c_char,len=ilen+1),pointer :: s
    call c_f_pointer(cp,s)
    fstr = s(1:ilen)
    nullify(s)
end block

end function c_ptr_to_f_string

The example function we are calling from Python is defined below. It is important to note that Fortran is expecting the cp argument to be passed by reference. That is why we use byref(cp) when we call it from Python.

subroutine send_json_to_fortran(cp) &
bind (c, name='send_json_to_fortran')

type(c_ptr),intent(in) :: cp
!! a `c_char_p` from python containing a JSON string.

type(json_file) :: json

call c_ptr_to_json(cp,json)

! do something with the data
! (in this case, just print it)
call json%print_file()

call json%destroy() ! free memory

end subroutine send_json_to_fortran

This prints the JSON structure (from Fortran) like so:

{
    "Generated in Python": true,
    "scalar": 1,
    "vector": [
        1,
        2,
        3
    ],
    "string": "hello"
}

Note that or send_json_to_fortran() routine has the BIND(c) attribute to make it callable from C (or in this case, Python). Once we have the JSON data on the Fortran side, we can also access the variables, query what's in the structure, get the values, etc., like we can any JSON data. If you want to modify the data on the Fortran side, keep reading.

Sending a JSON structure from Fortran to Python

The inverse operation (sending arbitrary data using JSON from Fortran back to Python) is also possible, but a little bit more complicated. First, in Fortran, we will define a container type that will hold a deferred-length (allocatable) string:

type,public :: container
character(len=:),allocatable :: str
end type container

The reason we need this is so that we can point to it with a pointer, which will be come apparent shortly. We then define a Fortran function to convert a JSON structure into a string, and then allocate a pointer to a container to hold the string:

subroutine json_to_c_ptr(json,cp,destroy)

implicit none

type(json_file),intent(inout) :: json !! JSON data
type(c_ptr) :: cp
!! a pointer to a container
!! containing the JSON data as a string
logical,intent(in) :: destroy
!! to also destroy the JSON file
!! (must be destroyed on the fortran
!! side somewhere to prevent memory leak)

character(len=:),allocatable :: str
type(container),pointer :: c

call json%print_to_string(str)
if (destroy) call json%destroy()
allocate(c)
c%str = str
cp = c_loc(c)

end subroutine json_to_c_ptr

We also will need a routine that is callable from Python to get the length of a string within a container:

function get_string_length(cp) result(ilen) &
bind (c, name='get_string_length')

type(c_ptr),intent(in) :: cp !! pointer to a container
integer(c_int) :: ilen !! the length of the string

type(container),pointer :: c

ilen = 0
if (c_associated(cp)) then
    call c_f_pointer (cp, c)
    if (allocated(c%str)) ilen = len(c%str)
end if

end function get_string_length

Now, we also need another routine callable from Python to populate a c_char_p string buffer:

subroutine populate_char_string(cp,string) &
bind(c,name='populate_char_string')

type(c_ptr),intent(in) :: cp !! pointer to a container
type(c_ptr),intent(inout) :: string
!! a preallocated string buffer that
!! the string will copied into

type(container),pointer :: c

if (c_associated(cp)) then
    call c_f_pointer (cp, c)
    if (allocated(c%str)) then
        call f_string_to_c_ptr(c%str, string)
    end if
end if

end subroutine populate_char_string

Finally, we also need a way to destroy the string that we are creating in Fortran. This has to be done in Fortran, since the container pointer was allocated in Fortran. Python will garbage collect its own variables, but it doesn't know how to garbage collect a Fortran pointer, so we'll have to call this routine from Python when we no longer need the variable.

subroutine destroy_string(cp) &
bind (c, name='destroy_string')

type(c_ptr),intent(inout) :: cp
!! pointer to a container

type(container),pointer :: c

if (c_associated(cp)) then
    call c_f_pointer (cp, c)
    if (allocated(c%str)) deallocate(c%str)
    deallocate(c)
end if

cp = c_null_ptr

end subroutine destroy_string

Now, say we have the following Fortran routine that generates some JSON data and sends it to Python:

subroutine get_json_from_fortran(cp) &
bind (c, name='test_send_json_to_python')

type(c_ptr) :: cp
!! pointer to a container containing a json string

type(json_file) :: json

! sample data:
call json%add('Generated in Fortran', .true.)
call json%add('scalar', 1)
call json%add('vector', [1,2,3])
call json%add('string', 'hello')
call json%add('string array', ['1','2','3'])

! convert it to a c_ptr (and destroy JSON structure)
call json_to_c_ptr(json,cp,destroy=.true.)

end subroutine get_json_from_fortran

On the Python side, we have:

# declare some more interfaces to the Fortran routines:
get_string_length = dll.get_string_length
get_string_length.argtypes = [POINTER(c_void_p)]
get_string_length.restype = c_int

populate_char_string = dll.populate_char_string
populate_char_string.argtypes = [POINTER(c_void_p),POINTER(c_char_p)]
populate_char_string.restype = None

destroy_string = dll.destroy_string
destroy_string.argtypes = [POINTER(c_void_p)]
destroy_string.restype = None

def fortran_str_to_python_dict(cp,destroy=True):

    # get the length of the string:
    length = c_int()
    length = get_string_length(byref(cp))

    # preallocate a string buffer of the correct size to hold it:
    s = c_char_p()
    s = cast(create_string_buffer(b' '.ljust(length)),c_char_p)

    # convert it to a normal python string:
    populate_char_string(byref(cp),byref(s))
    string = s.value.decode()

    if (destroy):
    # destroy it on the Fortran side:
    destroy_string(byref(cp))

    return json.loads(string)

The fortran_str_to_python_dict() routine is the main routine that converts the string that comes from Fortran into a Python dictionary. First we get the string length, then allocate a string buffer to hold it, we populate the buffer, and then parse the JSON. Note that we also have the option to destroy the Fortran string (which if you recall was allocated as a pointer in json_to_c_ptr() so it needs to be deallocated to prevent a memory leak).

Now we can use these routines from Python like so:

get_json_from_fortran = dll.get_json_from_fortran
get_json_from_fortran.argtypes = [POINTER(c_void_p)]
get_json_from_fortran.restype = None

cp = c_void_p()
get_json_from_fortran(byref(cp))
d = fortran_str_to_python_dict(cp,destroy=True)
print(json.dumps(d, indent=2))

Which prints the result:

{
    "Generated in Fortran": true,
    "scalar": 1,
    "vector": [
        1,
        2,
        3
    ],
    "string": "hello",
    "string array": [
        "1",
        "2",
        "3"
    ]
}

From Python to Fortran and Back

Finally, another use case would be to generate some data in Python, pass it to Fortran where it is modified, and then return it back to Python. In this case, we need to use the container approach, with a few additional routines. On the Python side:

c_ptr_to_container_c_ptr = dll.c_ptr_to_container_c_ptr
c_ptr_to_container_c_ptr.argtypes = [POINTER(c_char_p),POINTER(c_void_p)]
c_ptr_to_container_c_ptr.restype = None

def python_str_to_container(s):
    cp = python_str_to_fortran(s)
    ccp = c_void_p()
    c_ptr_to_container_c_ptr(byref(cp),byref(ccp))
    return ccp

def python_dict_to_container(d):
    return python_str_to_container(json.dumps(d))

And we also need a routine on the Fortran side to convert from the string buffer format to the container format:

subroutine c_ptr_to_container_c_ptr(cp,ccp) &
bind (c, name='c_ptr_to_container_c_ptr')

type(c_ptr),intent(in) :: cp
!! a `c_char_p` from python
type(c_ptr),intent(out) :: ccp
!! pointer to a container that contains the string

character(len=:),allocatable :: str
!! fortran version of the string
type(container),pointer :: c
!! container to hold the string

if (c_associated(cp)) then

    ! get the fortran string:
    str = c_ptr_to_f_string(cp)

    ! return a pointer to the container:
    allocate(c)
    c%str = str
    ccp = c_loc(c)

else
    ccp = c_null_ptr
end if

end subroutine c_ptr_to_container_c_ptr

Here is an example use case:

subroutine modify_json_in_fortran(cp) &
bind (c, name='modify_json_in_fortran')

type(c_ptr),intent(inout) :: cp !! a pointer to a container

type(json_file) :: json
type(container),pointer :: c

if (c_associated(cp)) then
    call c_f_pointer (cp, c)
    if (allocated(c%str)) then

        ! parse JSON:
        call json%load_from_string(c%str)

        !do something with the data:
        call json%print_file()
        call json%add('Added in Fortran', [9,10])

        ! convert it to a c_ptr (and destroy JSON structure)
        call json_to_c_ptr(json,cp,destroy=.true.)

    end if
end if

end subroutine modify_json_in_fortran

Which we can call from Python like so:

modify_json_in_fortran = dll.modify_json_in_fortran
modify_json_in_fortran.argtypes = [POINTER(c_void_p)]
modify_json_in_fortran.restype = None

cp = python_dict_to_container(a)
modify_json_in_fortran(byref(cp))
d = fortran_str_to_python_dict(cp,destroy=True)

print('')
print('modified by Fortran:')
print(json.dumps(d, indent=2))

Which prints:

{
    "Generated in Python": true,
    "scalar": 1,
    "vector": [
        1,
        2,
        3
    ],
    "string": "hello"
}

modified by Fortran:
{
    "Generated in Python": true,
    "scalar": 1,
    "vector": [
        1,
        2,
        3
    ],
    "string": "hello",
    "Added in Fortran": [
        9,
        10
    ]
}

Summary

Once all the interface routines are created, this technique is fairly easy to use and seems to work great. It is probably fine for a lot of use cases where extreme performance is not required. For example, it may be a good enough technique for exchanging data from a Python GUI and a core application written in Fortran. However, for very large data sets, this may not be very efficient since it does involve converting everything to and from strings.

See also

Updated on 7 June 2018: the code is now available on GitHub.

May 06, 2018

Low Lunar Orbit Maintenance

Altair_ascent_stage_docking_with_Orion

Altair's ascent stage docking with Orion [NASA]

The lunar gravity field is bumpy. Thus a spacecraft in a low lunar orbit for an extended period of time must perform periodic correction maneuvers or the orbit eventually becomes unstable and can crash into the Moon. When the Orion spacecraft was being designed during the late Constellation Program, the requirements of global surface access sortie missions, six-month polar outpost missions, anytime aborts, and autonomous operation of the vehicle while the crew was on the surface, required lots of analysis to figure out how to make the architecture work. And then the program was canceled.

The reference low-lunar orbit maintenance algorithm we used during Constellation was fairly simple. Orion initially inserted into a retrograde circular parking orbit with a 100 km altitude. During orbit propagation, if the altitude was reduced by a specified deadband altitude (say 10 km) the vehicle would perform a periapsis raise maneuver at the next apoapsis to bring the periapsis back up to 100 km. The apoapsis radius was not controlled, and indeed all the other orbital elements were also allowed to float. This would continue until the crew was ready to return (on the Altiar lander ascent stage), at which point Orion would circularize the orbit and perform a plane change to allow Altair to perform an in-plane ascent.

The \(\Delta v\) maneuver to perform at apoapsis of an orbit with a periapsis radius \(r_{pi}\) and apoapsis radius \(r_{ai}\), to target a periapsis radius \(r_{pf}\) can be computed by:

  • Initial apoapsis velocity: \(v_{ai} = \sqrt{ 2 \mu \frac{r_{pi}}{r_{ai} (r_{ai}+r_{pi})} }\)
  • Final apoapsis velocity: \(v_{af} = \sqrt{ 2 \mu \frac{r_{pf}}{r_{ai} (r_{ai}+r_{pf})} }\)
  • Maneuver magnitude: \(\Delta v = v_{af} - v_{ai}\)

We can implement this method using the Fortran Astrodynamics Toolkit fairly easily. We also need an integrator with an event finding capability (we can use DDEABM, a variable step size variable order Adams method). We propagate from the initial insertion state until the deadband altitude is violated, then propagate to the next apoapsis where the periapsis raise maneuver is performed, then repeat as necessary. The code implementing this scheme can be found on GitHub here. For the force model, we are using a 20 x 20 GRAIL geopotential model of the Moon (not available during Constellation, when we used the Lunar Prospector derived LP150Q model). A contour plot of the results for 100 km retrograde orbits with a 10 km altitude deadband and 10 days of propagation time are shown below. Each maneuver ends up being about 2 m/s for this example, and the bumpiness of the gravity field is evident. Some insertion orbits require no maneuvers, while others require about 18 m/s over the 10 day period.

lom

See also

Apr 15, 2018

New Coordinate Systems for Solar System Bodies

ann18010a

The International Astronomical Union (IAU) has published a new report with their latest recommendations for the cartographic coordinates and rotational elements of planetary bodies. The new guidelines, equations, and coefficients have been updated based on the latest data (including from interplanetary missions such as MESSENGER, Dawn, Rosetta, and Stardust-NExT).

One difference from the previous reports is that the low-fidelity Earth and Moon frames have been removed "in order to avoid confusion", since high-accuracy Earth and Moon frames are available. I mentioned the IAU_EARTH and IAU_MOON frames in a previous post. They are available in SPICE and the Fortran Astrodynamics Toolkit. I think it's unfortunate that they were removed. It's true that these are not ultra-precise frames (see, for example, this tutorial from JPL). However, low-to-medium fidelity frames are still very useful for many applications, including spacecraft trajectory optimization. A frame that is good enough, is very fast to compute, is continuous and differentiable for all time, and can be implemented in a few lines of code is often preferable to having to deal with thousands of lines of SOFA or SPICE code, leap second confusion, binary kernel files, etc. that are required to get the high-accuracy frames to work.

See also

Apr 04, 2018

Orbital Mechanics with Fortran

crtbp

Many of NASA's workhorse spacecraft trajectory optimization tools are written in Fortran or have Fortran components (e.g., OTIS, MALTO, Mystic, and Copernicus). None of it is really publicly-available however. There are some freely-available Fortran codes around for general-purpose orbital mechanics applications. Including:

You can also find some old Fortran gems in various papers on the NASA Technical Reports Server.

See also

  1. J. Williams, R. D. Falck, and I. B. Beekman. "Application of Modern Fortran to Spacecraft Trajectory Design and Optimization", 2018 Space Flight Mechanics Meeting, AIAA SciTech Forum, (AIAA 2018-1451)

Mar 10, 2018

OpenCoarrays 2.0.0 Released

PGAS

OpenCoarrays 2.0.0 has been released. According to the release notes, this release adds support for:

  • Fortran 2018 teams (note that this requires the use of GFortran 8.x, currently under development)
  • Fortran 2018 QUIET specifier in STOP and ERROR STOP statements
  • Fortran 2008 allocatable & pointer components of derived-type coarrays

Since Fortran 2008 (with the inclusion of coarrays in the standard), Fortran has been a PGAS programming language. Fortran 2018 adds additional enhancements for coarrays (including the teams feature, and support for detection of failed images). OpenCoarrays is an open source (permissive BSD license) library that enables coarray support in the GFortran compiler.

See also

Feb 17, 2018

Command Line Arguments

prompt

A standard way to access a program's command line arguments wasn't added to Fortran until 2003. Command line arguments were one of those things (like environment variables, the file system, and the operating system shell) that the Fortran standards committee pretended didn't exist or wasn't important. Of course, that was and is ridiculous. So for decades, all Fortran compilers provided a non-standard way to get the command line arguments (for example, GETARG in gfortran). Currently, the standards committee still pretends that the file system doesn't exist (for example, there's no standard way to change an application's working directory, or list the files in a directory, or even to check for the existence of directories). But I digress.

The two routines you need to know are COMMAND_ARGUMENT_COUNT and GET_COMMAND_ARGUMENT. Frequently in examples on the internet (even the gfortran example) you will see people just declare a really long string and use GET_COMMAND_ARGUMENT to get the command line argument and hope for the best. That's terrible of course, since this routine does have the ability to return the actual length of the argument, so you should use that to allocate a string with the correct size to hold it.

Consider the following example:

module command_line_module

    implicit none

    private

    type,public :: string
        !! a variable-length string type
        character(len=:),allocatable :: str
    end type

    public :: get_command_arguments

contains

subroutine get_command_arguments(nargs, args)

    !! return all the command line arguments.
    !! (each one is the correct size)

    implicit none

    integer,intent(out) :: nargs !! number of arguments
    type(string),dimension(:),allocatable,intent(out) :: args

    integer :: i, ilen

    nargs = command_argument_count()
    allocate(args(0:nargs))

    do i = 0, nargs
        call get_command_argument(i,length=ilen)
        allocate(character(len=ilen) :: args(i)%str)
        call get_command_argument(i,value=args(i)%str)
    end do

end subroutine get_command_arguments

end module command_line_module

Here we have to define our own string type (string), so we can create an array of deferred-length strings (since a useful varying-length string type is another thing Fortran pretends is not important). See StringiFor for a full-featured variable-string type, but for our purposes here the simple one is fine. The procedure in the get_command_arguments() routine is to check for the size of each argument, allocate a string large enough to hold it, and then retrieve the value. It returns an array of string variables containing the command line arguments. With a wrapper module like this, getting command line arguments is effortless, and there is no reason anymore to use the old non-standard routines.

Feb 10, 2018

Another One Bites the Dust: GRAM Atmosphere Model

It looks like the Fortran community has lost another venerable NASA Fortran library to C++. I recently noticed that the latest release of NASA's Global Reference Atmosphere Model (GRAM) is now C++. From the release link:

Earth-GRAM 2016 is now available as an open-source C++ computer code that can run on a variety of platforms including PCs and UNIX stations. The software provides a model that offers values for atmospheric parameters such as density, temperature, winds, and constituents for any month and at any altitude and location within the Earth's atmosphere. Earth-GRAM 2010 is available in FORTRAN. Similar computer models of Mars, Venus, Titan, and Neptune are available for diverse mission applications through NASA’s Space Environments and Effects (SEE) Program, a partnership with industry, academia, and other government agencies that seeks to develop more reliable, more effective spacecraft.

globe_west

The Earth-GRAM model, originally written in Fortran 77, is produced by MSFC and has been around for a while, with many revisions (i.e., GRAM-86, GRAM-88, GRAM-90, GRAM-95, GRAM-99, GRAM 2007, and GRAM 2010). In the 2010 revision, they converted the old Fortran 77 code to Fortran 90, which I thought was a step in the right direction (although maybe 20 years too late). It seems like they could have used Fortran 2003 and provided a C interface using the standard Fortran/C Interoperability feature, which would have enabled calling it from a variety of languages (including C++).

I remain utterly unconvinced that C++ is the programming language that engineers and scientists should be using. Are CS people taking over programming jobs once done by engineers? Is this a failure of university engineering departments? Although it seems like even CS people don't want to be using C++ anymore. Every time I turn around a new programming language (Java, Rust, Go, D, C#, Swift, ...) is created by CS people fed up with the C++ dumpster fire. So I don't know who is telling engineers that C++ is the future. And the perception that Fortran is an obsolete programming language is very strong, even within the engineering/technical disciplines that it was designed for. Sometimes this is due to sheer ignorance, but sometimes it is due to ancient Fortran 77 legacy code that people have to deal with. I sympathize with them, Fortran 77 code is terrible and needs to go away, but it should be modernized, not thrown out and replaced with C++. Think of the children! (and the memory leaks!)

Parting words of wisdom

Security tips when programming in C (2017 edition):

1) Stop typing

2) Delete what you've already typed

— ryan huber (@ryanhuber) June 21, 2017

See also

← Previous Next → Page 4 of 14