Programming challenge (19th century writers)

Personally, I’d also be happy with a subroutine:

call order_inplace(table,indx)

I can imagine delegating both the partitioning and reordering to a derived type which accepts a predicate function, or a mask expression. I’m just not sure how to make it sufficiently generic?

Could a generic reorder be achieved with a class(*) input argument? After all, at this point we don’t care about the contents of the table items but just that they are moved to the right place.

Edit: with respect to partioning

The partitioned subroutine, could use the true_pos function of @beliavsky given here.

This was rather clever since the age array is temporary anyways. It does require two searches though.

I think it can be made neater with associate

  subroutine find_second_oldest( table )
    ! Determine the second oldest person as quickly as possible ('2')
    ! Note that you are allowed to change the order of persons.
    type(person), intent(inout) :: table(:)
    integer :: age(size(table))
    age = table%death_year - table%birth_year  ! we could also introduce a function
    age(maxloc(age,dim=1)) = -1
    i = maxloc(age,dim=1)
    associate( p => table(i))
      print *, "second oldest age = ", age(i)
      write(*,prnt) p%firstname, p%lastname, p%birth_year, p%death_year
    end associate

  end subroutine

Maybe a function like nmaxloc would make a nice addition (as well as nminloc, nminval, nmaxval) to stdlib. I’m not sure if these would better be done with partial sorting though, instead of n searches?

Even nicer syntax could be achieved with an (impure?) function called nth_element:

associate( p => table( nth_element(2,age(table),back=.true.) ) )
  print *, "second oldest age = ", age(p)
end associate

elemental integer function age(p)
  type(person), intent(in) :: p
  age = p%death_year - p%birth_year
end function

Same goes here, an associate helps reduce the verbosity:

associate(mask => table%birth_year.gt.1900 .and. table%birth_year.lt.2000)
    table=[ pack( table, mask), &
          & pack( table, .not. mask ))]
   ! alternatively
   ! call stable_partition_inplace(table,mask)
end associate

Maybe a new array transformational function?

partition = [pack(x,mask),pack(x,.not. mask)]  ! incomplete example

Some intrinsic functions for the operations done in STL would really enhance the Fortran programming experience.

I made a few of the changes in the original, as you can see the diff by clicking on the history (instead of a second post). Was looking for a way to let anyone edit; not sure if that is the default or not. A lot of the solutions are just a line or two and as mentioned could just be done in the SELECT but it seemed more the intent to stick to the original structure so I did not do that (but it would make the code very short).

Some interesting thoughts on better Fortran or stdlib capabilities. It got deeper and more interesting than I anticipated (which is great!). Tempted to stick this in a github repo and use the commit log to capture some of the changes and start a collection of examples but not feeling that ambitious, and already have about forty repos that need finished already :wink:

1 Like

The lecturer in the course is one of the best C++ educators/consultants, so I think he deserves credit for the format of the exercise. It gave me lots of thoughts on how to approach Fortran, but I still have to reflect on them for a while.

Although you could achieve many things with one- or two-liners in the main program, that is not the goal. The idea is to follow Separation of Concerns (SoC). By handling the operations in sub-programs, one can immediately see the main program structure and intent with the repeated do loop and operations, while all the “nasty” details are invisible. It gives the source code a nice sense of cleanliness and symmetry, making it easy to grasp and extend.

In fact the original exercise handled 8 cases and I reduced them for this thread. I was excited to see you added the “bonus” cases back. :slight_smile:

Well, I have not achieved complete genericity, but with a bit of help, I do have a subroutine that can reverse an array in place without taking particular care about the data type:

! reverse_array.f90 --
!     Attempt to create a generic routine that reverses an array - in place
!
module reverse_array
    implicit none

contains
subroutine reverse( array, assign )
    class(*), dimension(:), intent(inout), pointer     :: array
    interface
        subroutine assign( a, b )
            class(*), intent(out) :: a
            class(*), intent(in)  :: b
        end subroutine assign
    end interface

    class(*), allocatable                              :: tmp
    integer                                            :: i, j

    allocate( tmp, mold=array(1) )

    do i = 1,size(array)/2
        j        = size(array) + 1 - i
        call assign( tmp, array(i) )
        call assign( array(i), array(j) )
        call assign( array(j), tmp )
    enddo
end subroutine reverse
end module reverse_array

! test --
!     Test this implementation
!
program test_reverse_array
    use reverse_array

    type :: person
        ! Each person has a firstname, lastname, birth_year, and death year.
        character(len=30) :: firstname
        character(len=30) :: lastname
        integer           :: birth_year
        integer           :: death_year
    end type

    type(person), allocatable, dimension(:), target :: table
    class(*), dimension(:), pointer                 :: ptable

    table = [ &
        person( "Miguel",  "de Cervantes", 1547, 1616 ), &
        person( "Charles", "Dickens",      1812, 1870 ), &
        person( "George",  "Orwell",       1903, 1950 ), &
        person( "Ernest",  "Hemingway",    1899, 1961 ), &
        person( "William", "Shakespeare",  1564, 1616 ), &
        person( "Leo",     "Tolstoy",      1828, 1910 ), &
        person( "Scott",   "Fitzgerald",   1896, 1940 ), &
        person( "Franz",   "Kafka",        1883, 1924 ), &
        person( "Mark",    "Twain",        1835, 1910 ), &
        person( "Roald",   "Dahl",         1916, 1990 ), &
        person( "James",   "Joyce",        1882, 1941 ) ]

    ptable => table
    call reverse( ptable, assign_person )

    write(*,'(2a,2i6)') table

contains
subroutine assign_person( a, b )
    class(*), intent(out) :: a
    class(*), intent(in)  :: b

    select type( a )
        type is (person)
            select type( b )
                type is (person)
                    a = b
            end select
    end select
end subroutine assign_person
end program test_reverse_array

Out of laziness, I stole, eh, copied, the person data type and the table of persons from the above example.

1 Like

If the code is modified to also work with a date derived type, it looks as shown below. The subroutines assign_person and assign_date are easy to write but also seem like boilerplate. I don’t really understand OOP in any language, but could Fortran be extended so that there was a default assign subroutine for two arguments that are objects of the same type? Or is there already a way to write a generic assign subroutine?

code
! reverse_array.f90 --
!     Attempt to create a generic routine that reverses an array - in place
!
module reverse_array
    implicit none

contains
subroutine reverse( array, assign )
    class(*), dimension(:), intent(inout), pointer     :: array
    interface
        subroutine assign( a, b )
            class(*), intent(out) :: a
            class(*), intent(in)  :: b
        end subroutine assign
    end interface

    class(*), allocatable                              :: tmp
    integer                                            :: i, j

    allocate( tmp, mold=array(1) )

    do i = 1,size(array)/2
        j        = size(array) + 1 - i
        call assign( tmp, array(i) )
        call assign( array(i), array(j) )
        call assign( array(j), tmp )
    enddo
end subroutine reverse
end module reverse_array

! test --
!     Test this implementation
!
program test_reverse_array
    use reverse_array, only: reverse

    type :: person
        ! Each person has a firstname, lastname, birth_year, and death year.
        character(len=30) :: firstname
        character(len=30) :: lastname
        integer           :: birth_year
        integer           :: death_year
    end type person
!
    type :: date
        integer :: year, month, day
    end type date

    type(person), allocatable, target :: table(:)
    type(date)  , allocatable, target :: dates(:)
    class(*)    , pointer             :: pvec(:)
    dates = [date(2021,1,1),date(2020,1,1),date(2019,1,1)]

    table = [ &
        person( "Miguel",  "de Cervantes", 1547, 1616 ), &
        person( "Charles", "Dickens",      1812, 1870 ), &
        person( "George",  "Orwell",       1903, 1950 )]

    pvec => table
    call reverse( pvec, assign_person )

    write(*,'(2a,2i6)') table
    pvec => dates
    write (*,"(' initial:',*(1x,i0))") dates
    call reverse( pvec, assign_date )
    write (*,"('reversed:',*(1x,i0))") dates
contains
subroutine assign_person( a, b )
    class(*), intent(out) :: a
    class(*), intent(in)  :: b
    select type( a )
        type is (person)
            select type( b )
                type is (person)
                    a = b
            end select
    end select
end subroutine assign_person
!
subroutine assign_date(a,b)
class(*), intent(out) :: a
class(*), intent(in)  :: b
select type (a)
   type is (date)
      select type(b)
         type is (date)
            a = b
      end select
end select
end subroutine assign_date
end program test_reverse_array
ouptut
George                        Orwell                          1903  1950
Charles                       Dickens                         1812  1870
Miguel                        de Cervantes                    1547  1616
 initial: 2021 1 1 2020 1 1 2019 1 1
reversed: 2019 1 1 2020 1 1 2021 1 1

This sure looks like a template. Would be pretty nice to figure out which “default” templates are needed to work with a family of generic reordering routines.

On a different note, would the following also work?

subroutine assign_person( a, b )
    class(*), intent(out) :: a
    class(*), intent(in)  :: b

    if (same_type_as(a,b))  a = b

end subroutine assign_person

Not at present, but maybe you are suggesting a language extension. Gfortran says

   96 |     if (same_type_as(a,b))  a = b
      |                            1
Error: Nonallocatable variable must not be polymorphic in intrinsic assignment at (1) - check that there is a matching specific subroutine for '=' operator

Yes, that sort of messages is the reason I pass a pointer to the array, instead of the array itself.

Perhaps transfer() can be persuaded to help out?

I guess the same_type_as can’t work in this case, because the concrete type is still unknown, and only becomes known in the select type block…

Thank you all for the very interesting discussion. I thought of dropping a quick message noting that some of what you have shown here could probably be used as the skeleton for methods of a stdlib data frame. Some simple implementations exist already.

1 Like

The conversation got interesting along other avenues so I was not going
to post these changes and interrupt things, but I was making a version
to play with some of the ideas expounded here and along the way decided
to see if I could make procedures conforming to the rules that I had
used from fpm(1) dependencies.

I made a one-line sort(3f) (maybe a known method, but I had not seen
it and was rather pleased with it but have not tested it much yet :>
). What really surprised me, and is related to these discussions is that
when I used it in the form

   table( sort(table%lastname) )=table

it worked fine with gfortran(1) and ifort(1), but when I used

   table=table( sort(table%lastname)

bizarre errors occurred, including ifort(1) apparently executing it twice.
I don’t want to distract from the discussions here but if anyone has time
to see if that sort(3f) is doing something illegal let me know. I don’t see
it so far.

git clone https://github.com/urbanjost/chal.git
cd chal
fpm run chal
# or you can just compile chal/app/main.f90 (no dependencies at the moment)
# as it is just one file; but fpm(1) rules!

Thanks @epagone for the interesting perspective. I didn’t notice it before, but indeed the table and the operations that can be performed are reminiscent of a spreadsheet or dataframe.

I must admit, I don’t really see how a Fortran dataframe could be superior to what is available in R, Pandas, or Excel. These tools allow interactive usage, i.e. you display some data, evaluate, modify your commands, and reprint. Essentially you work with a REPL (in case of Excel a GUI).

Currently, a Fortran program is compiled, meaning unless you write a full blown expression parser, you can only support a very limited set of operations. Sure it can be nice to take a peak at some data on the fly, but I think the other available tools are far superior to what can be done in Fortran. Undertaking a project like this feels to me like re-inventing the wheel. But I can imagine with LFortran providing a REPL, a dataframe will begin to make sense.

Does the behavior match in both cases, if you write the indexes to an array first, instead of providing them directly as an expression?

I have many Fortran programs that work with daily times series of stock prices or returns. Instead of passing around separately dates(:), returns(:,:), stock_symbols(:) to procedures and checking in each procedure whether they are dimensioned consistently I use a derived type

type :: date_frame_full
   character (len=1000)                 :: title = ""
   type(date_mdy)         , allocatable :: dates(:)  ! (nobs)
   character (len=len_sym), allocatable :: sym(:)    ! (nvar)
   real(kind=dp)          , allocatable :: xx(:,:)   ! (nobs,nvar)
end type date_frame_full

for which I have defined

many procedures
public :: date_frame_full,set_date_frame,read_date_frame,num_var,num_obs, &
          num_obs_total,first_date,last_date,compute_returns,write_date_frame,bad_real, &
          date_frame_of_file,first_year,resample_period, &
          weighted_sum,allocate_from_df,random_df,lag,operator(-),operator(+),mult,select_var, &
          conform,write_date_frame_stream,alloc_df,read_date_frame_stream, &
          serialize_date_frame,select_obs,shape_df,reading_time, &
          print_shapes,df_product,df_rows_sumproduct,add_df,df_skip_initial_obs,display, &
          display_tail,display_head_tail,compute_returns_with_dividends,dates_pos,var_pos,project_dates, &
          project_var,project_dates_var,obs_since,cumul_product,operator(/), &
          operator(*),project_date_frame,date_range_string,demean_by_var,first_few_var,pos_date,dates_df_vec, &
          subset,add_var_vec,add_dividend_adjusted_prices,add_returns,fill_data,subset_obs_bounds,operator(**), &
          abs,subset_sym

I have found it useful and could clean it up and publish it if there is interest. One thing I have wondered is how to create a similar derived type

type :: date_frame_seconds ! time measured to integer seconds
   character (len=1000)                 :: title = ""
   integer                , allocatable :: times(:)  ! (nobs) -- seconds after midnight
   character (len=len_sym), allocatable :: sym(:)    ! (nvar)
   real(kind=dp)          , allocatable :: xx(:,:)   ! (nobs,nvar)
end type date_frame_seconds

without excessive code duplication. Such derived types could also be used for time series in climatology and other domains.

1 Like

I reported the definite bugs; but the remaining one is just me, I think. Apparently something I thought was equivalent is not, as I get the same error with three compilers. I will figure it out when I have time. The exercise has been fun and productive, but I am apparently having a brain burp. I thought both of these sorts would work. I know it will be something simple but I am not seeing it:

program main
! bug: ifort calls sorti(1) twice for each invocation
implicit none
integer,parameter            :: isz=100
real                         :: rr(isz)
integer                      :: ii(isz), jj(isz)
integer                      :: i
character(len=*),parameter :: x='(*(g0,1x))'
   !write(*,*)'HUGE=',huge(0) !  HUGE=  2147483647
   call doboth( [100,100,500,300,100,400,500,200] )
   call doboth( [0,0,0,0,0] )
   call doboth( [100,500,300,400,200] )
   call doboth( [100,500,-300,0,400,200] )
   call doboth( [600,300,200,100,0,-100,-200,-300] )
   write(*,x)repeat('=',80)
   write(*,*)'initializing array with ',isz,' random numbers'
   CALL RANDOM_NUMBER(RR)
   jj=rr*450000.0
   ! use the index array to actually move the input array into a sorted order
   ii=jj
   ii=ii(sorti(ii))
   write(*,*)'checking if values are sorted(3f)'
   do i=1,isz-1
      if(ii(i).gt.ii(i+1))then
         write(*,*)'Error in sorting reals small to large ',i,ii(i),ii(i+1)
      endif
   enddo
   write(*,*)'test of sorti(3f) complete'

   ii=jj
   ii(sorti(ii))=ii
   write(*,*)'checking if values are sorted(3f)'
   do i=1,isz-1
      if(ii(i).gt.ii(i+1))then
         write(*,*)'Error in sorting reals small to large ',i,ii(i),ii(i+1)
      endif
   enddo
   write(*,*)'test of sorti(3f) complete'
contains
subroutine doboth(ints)
integer :: ints(:)
integer,allocatable :: iii(:)
   iii=ints
   write(*,x)repeat('=',80)
   write(*,x)'INPUT',iii
   iii(sorti(iii))=iii
   write(*,x)'LHS  ',iii
   iii=iii(sorti(iii))
   write(*,x)'RHS  ',iii
end subroutine doboth
function sorti(ints) result(counts)
integer :: ints(:)
integer :: counts(size(ints)), i
integer, save :: calls=0
   ! WARNING: this sort is slow as mollasses, but if fun to write in one line
   ! relative sort first            take account of duplicates
   counts=[(count(ints(i) > ints) + count(ints(i) == ints(:i)), i=1,size(ints) )]
   calls=calls+1
   write(*,x)'CALLS',CALLS,'VALUES',counts
end function sorti
end program

I added a few things that push modern features and even though I was totally content with using fpm(1) dependencies I made a stand-alone version in app/main.f90 if any non-fpm users want to start with my version (which does not use the logic above with the bug) in

if anyone else wants to play. Ignore the bug.f90 and blocky.f90 files. So good for now (unless someone wants to explain why my first sort in above does not work :>).

There are some interesting ideas discussed above that would be interesting to implement. I am curious to compare the app/main.f90 version to some of the C++ versions if any classmates are brazen enough to post in hostile terroritory :smile:

So several bug reports, and the M_random and chal github repositories and some fun posting later I am taking a break from worrying about it for now, but I learned I did not know something I thought I did, which is always useful. I originally just thought it was a good puzzle to play with while I drank a cup of coffee.

PS: concerning some of the comments regarding something like a spreadsheet or shell in Fortran mentioned above:

I have used a program for years (as have several other hundred people) that is a fortran shell with a full expression parser, graphics, steam tables, integration, differentiation, several command styles including a purely functional one (even IF is a function, but it can execute named blocks of code as well as expressions) that I still prefer over gnuplot, R, MATLAB, and so on. Even today a big advantage (it has not changed radically since the early 90’s) is that you can load Fortran functions into it directly without creating bindings; and the people that use it have a LOT of those, and the data formats (binary and text) that it can read are very simple for Fortran programs to write, so Fortran is a totally adequate language for that type of stuff, it’s just that most places have not maintained or developed them. Some other old tools are still around and have an impressive amount of functionality, I also used a nearly complete FORTRAN77 interpreter interactively for years, so I am not such a hater of implicit typing as some are. Interactively, it is quite handy. Unfortunately, that is gone (I think). I still don’t know anything else where I can say to take the output of two runs and produce deltas on temperature versus pressure tables, and anywhere where the delta is > a threshhold make automatically labeled plots showing the curves to the left axis and the delta to right axis with a single command, and play back everything including mouse selections from a journal file (which can also be used to verify exactly what operations were performed on the data) in any of those other tools (and pedigree is critical to me) so I don’t think I am going to switch anytime soon.

Sounds interesting – I assume it is proprietary. Can parts of it be recreated in LFortran? Even if the code or an executable cannot be released, releasing the documentation could be beneficial. Sometimes I look at the documentation Matlab, NAG, IMSL, and other commercial software for algorithms and references.

It is interesting, at least to me, but I and others have tried several times and been repeatedly told no. The company involved is not in it’s prime and there has not been a new release for nearly six years now, which is sad to me. The syntax is Unix-like (actually, AEGIS-like, but that is a long story – we had never used Unix when it was started) so a lot of people familiar with Unix/GNU-Linux find it easy to use, and like other shells any command not built in is passed to the OS so it feels like you are in a typical shell when using it. You can use “cd” and “ls” and so on, for example (well, cd is actually a built-in but a look-alike of other shells). It lets you build custom commands somewhat like bash functions but you define a command prototype so you get automatic command line parsing. The journaling is particularly handy as it includes all non-graphic output with screen messages recorded as comments (it ALWAYS journals) so you can replay anything you did. The functional mode is something that is an acquired taste so not everyone uses it, but
you can do interesting things that are somewhat like SQL in functionality with it. It is usually used for extracting data from tables, but it can be powerful.Something like

plot -f 3 -if eq(max(ge(temp,300) ge(pressure,400),primary(“HG1*”),0)

is one of the uglier but powerful commands you can use. That is not typical. Sort of like RPN.
Some people love it, others hate it and use a series of simpler commands instead. A main initial design feature was to make it feel like the sh(1) shell and to use Fortran syntax for expressions, which a lot of the targeted users were familiar with so you could start using it quickly. Never did add pipes to it though, which would have been a great addition. But you use files like tables so you don’t appear to be using files for data much except to load them and save them.

Because the primary file types are self-describing labeling and units are known so a lot of the users use just a few commands like plot, math, attach, save, and read; but others made some pretty impressive programs with it. I remember some that took projects projected to take several years and automated them and ran them in a few days.

Some of the code used for the prototype are publicly available but really dated. I could probably paste M_kracken, M_graph, M_calculator together and re-create the first production version (circa 1990-ish!) but that would not be very satisfying.

That was a fun code to work on, unfortunately since the Tsunami (long story why that is related) I don’t think anyone has worked on it except to port it to new machines.

Just one screen from the TOC for the help (there are about ten pages of main topics, some like steam(1) and matrix(1) go to many others:

paranoia (1ush)      - [USH-D:CMD] test floating point operations
parcel (1ush)        - [USH-D:CMD] create a new command
relate (1ush)        - [USH-D:CMD] unit conversion
setdash (1ush)       - [USH-D:CMD] create dash-code patterns for drawing lines
setmark (1ush)       - [USH-D:CMD] create or change geometric markers 1 thru 20
show (1ush)          - [USH-D:CMD] display command keywords and values
system (1ush)        - [USH-D:CMD] turn access to system commands on and off
units (1ush)         - [USH-D:CMD] control or display text associated with unit codes
wipe (1ush)          - [USH-D:CMD] empty the pseudo file (file 0)
zoom (1ush)          - [USH-D:CMD] select X11 keyboard zoom mode for plot(1) command
ansi (3ush)          - [USH-E:CALC] basic Fortran-intrinsic calculator functions
bessel (3ush)        - [USH-E:CALC] Bessel functions in the calculator
c (3ush)             - [USH-E:CALC] functions for mixing curves from multiple files
Calculator (3ush)    - [USH-E:CALC] using the calculator in USH
clcon (3ush)         - [USH-E:CALC] color model conversion calculator function
logical (3ush)       - [USH-E:CALC] NUMERIC LOGICAL FUNCTIONS
prompting_functions (3ush) - [USH-E:CALC] single numeric OR string prompting functions
steam (3ush)         - [USH-E:CALC] ASME steam functions
string (3ush)        - [USH-E:CALC] string functions
time (3ush)          - [USH-E:CALC] Time functions in jucalc(3f) expressions
wif97 (3ush)         - [USH-E:CALC] IF-97 Steam Routines
X11_functions (3ush) - [USH-E:CALC] X11 mouse functions
variables (7ush)     - [USH-E:CALC] miscellaneous functions and variables
convert (3ush)       - [USH-F:OPERATOR] A linear conversion operator for math(1)
decimate (3ush)      - [USH-F:OPERATOR] extract a portion of a curve
dif (3ush)           - [USH-F:OPERATOR] derivative of a curve
intg (3ush)          - [USH-F:OPERATOR] integrate curves
lfit (3ush)          - [USH-F:OPERATOR] linear fit of curves
lowess (3ush)        - [USH-F:OPERATOR] weighted regression smoothing operator for math(1)
map (3ush)           - [USH-F:OPERATOR] map curves to new X values
mask (3ush)          - [USH-F:OPERATOR] select curve points with a mask
operators (3ush)     - [USH-F:OPERATOR] OPERATORS for get(1) and math(1) (calculus, fitting, data reduction, ..
poly (3ush)          - [USH-F:OPERATOR] Polynomial curve fit of order 1 <= n <= 10.
page (1ush)          - [USH-G:LOWLEVEL] initializing the display for low-level graphics mode
parea (1ush)         - [USH-G:LOWLEVEL] creating subplot areas in low-level graphics mode
2 Likes

Actually, the pointer is not necessary. It was a relic from an earlier attempt. But unfortunately, using transfer() is not an option, as it runs into the same problem as attempting the switch directly.