I made a calculator/Interpreter in Fortran that works with scalars and 1D arrays of reals. It has a bit of the functionality of NumPy or R. You can define variables and evaluate expressions involving them. There are no loops or conditionals. There are whole array operations, array sections, vector subscripts, uniform and normal RNG, descriptive statistics, sorting, relational operators, and functions for correlation and covariance. ChatGPT o3 is helpful in writing a parser. One can use it to look up values for special functions, for example
To simulate the probability that a standard normal deviate exceeds 2 one writes
x = rnorm(10^6)
size mean sd min max first last
1000000 -0.0002 1.0001 -4.6193 4.9802 -0.8247 -0.0022
mean(x > 2)
0.022749
The last line works since logical variables are treated as 1s and 0s (not standard Fortran).
I don’t know how easily the calculator can be extended to handle 2D arrays. (M_matrix of @urbanjost does.) By toggling a parameter mutable you can allow or disallow variables to be changed once set. It’s mostly a fun project, but it is nice to quickly define functions for an interpreted language in a few lines of Fortran. It compiles with either gfortran or ifx.
I want to keep ^ for exponentiation and am not sure that ** should also do the same thing. Maybe Fortran used ** historically because ^ was not universally available. In the main program there is now code
line = replace(line, "**", "^")
call eval_print(trim(line))
so that ** is replaced by ^ before the code is processed. Now you get
Very interesting! Incidentally, I saw that you coded the function cumsum in one of the project modules
function cumsum(x) result(y)
! return the cumulative sum of x
real(kind=dp), intent(in) :: x(:)
real(kind=dp) :: y(size(x))
integer :: i, n
n = size(x)
if (n < 1) return
y(1) = x(1)
do i=2,n
y(i) = y(i-1) + x(i)
end do
end function cumsum
I copied it and tested it, but for large arrays it gives a stack overflow. I prefer not using the compilation flag heap-arrays0 (it might degrade performance?) so I made the function result as allocatable and it works without errors.
I wonder if it is a good practice to declare function results as allocatable to avoid stack overflow. Are there better alternatives (besides heap arrays flag)? Generally, I always try to avoid “automatic arrays” (not sure if this term is correct) in functions and subroutines because sooner or later my code crash
function cumsum(x) result(y)
! return the cumulative sum of x
real(kind=dp), intent(in) :: x(:)
real(kind=dp), allocatable :: y(:)
integer :: i, n
n = size(x)
allocate(y(n))
if (n < 1) return
y(1) = x(1)
do i=2,n
y(i) = y(i-1) + x(i)
end do
end function cumsum
The programmer should have a way to tell the compiler to first try stack allocation, and then if that fails to fall back to heap allocation for automatic arrays. That way, the program only does the extra overhead when necessary. However, I don’t know of any compilers that work that way. I don’t know why.
Nice!
The advantage of ^ over ** is that it keeps the parser simple and clean since you can proceed one character at a time. ** requires adding some extra logic to check the i+1 character and take care of exceptions (e.g. make sure the i+1 index does not exceed the string length, etc.).
If you want to extend the calculator to matrices you can also have a look at what @fluidnumerics_joe did in feq-parse. He handles up to 4 dimensions. If all the functions are elemental this is ‘just’ about creating a dedicated stack for 2D arrays. There is also some support to add user defined functions that could be fun to play with.
I changed cumsum to use an allocatable result. I also added the merge, pack, and count intrinsic functions, added print_compiler_version and print_compiler_info functions that call the compiler_version and compiler_option intrinsics, added a time_code toggle in the main program that times each line of code, and fixed some bugs.
You might have meant fortran compiler and not compilers in general, but at least one Ada compiler does something similar with a heap-allocated “secondary stack”: https://nytpu.com/gemlog/2024-12-27-2. And I previously read this on Wikipedia under the heading “System Interface”:
A safer version of alloca called _malloca , which allocates on the heap if the allocation size is too large, and reports stack overflow errors, exists on Microsoft Windows. It requires the use of _freea . gnulib provides an equivalent interface, albeit instead of throwing an SEH exception on overflow, it delegates to malloc when an overlarge size is detected.
Seems like something like that could be used to implement what you suggested (though maybe only on GNU and Windows).
I think most compilers allocate on a heap using some heap allocator, but I also like the secondary stack (per thread) better. I think you only need to manage it in functions that can allocate on it. So it should be more performing than a heap.
That approach should work also for fortran. The main advantage of stack storage is its simplicity, which results in its efficiency, over more general heap storage. It should make little difference in efficiency if it is a stack managed by the OS or a stack managed by the compiler. Of course, heap allocations cannot be eliminated entirely, as push/pop stack allocation only satisfies a subset of memory allocation requirements within the language (candidates are things like local automatic arrays, local allocatable arrays, allocatable function results, and temporary arrays needed during expression evaluation).
I added plotting with gnuplot (but if you substitute your own module plot_mod with the same subroutine signatures, other plotting tools can be used). An example is
There are one-line if statements and do loops with cycle and exit, so one can write
do i=1,10
x = rnorm()
if (x < 0) cycle
[x sqrt(x)]
end do
Sometimes you just want to repeat a line several times and do not need a loop variable. The interpreter treats *n at the beginning of a line as an instruction to execute the line n times, so
a = 3
x = 1
*5 x = (x + a/x)/2 ! Newton's method for square root
gives
> 3
> 1
> 2
1.750000
1.732143
1.732051
1.732051
and
*3 rnorm(10^6) ! simulate 3 sets of a million normal deviates
gives
size mean sd skew kurt min max first last
1000000 -0.0003 0.9995 -0.0024 0.0009 -5.9611 5.0877 0.8383 1.8332
size mean sd skew kurt min max first last
1000000 0.0008 1.0000 0.0006 0.0018 -4.6125 4.7730 0.3212 0.6562
size mean sd skew kurt min max first last
1000000 -0.0008 0.9993 -0.0005 -0.0017 -4.8507 5.3175 -0.7557 1.4073
