4 ms·
I really don't understand why you're arguing with my explanation. I'm not attacking fortran, I'm explaining why it will feel old to people: because it has thing
by BlackFly 5y ago
I really don't understand why you're arguing with my explanation. I'm not attacking fortran, I'm explaining why it will feel old to people: because it has things in it that people moved on from. Heck, I'd even agree that it is a sort of bad attitude to automatically discount something because of such legacy. But there are clear historical artifacts in the language (like most languages today).
> A very clever way of applying your algo to a submatrix from the days of sequence association and no array expressions.
This isn't why you specify the leading dimension, it is just a useful feature of the leftover artifact of the lack of dynamic allocation. If the algorithms were written today, the array slicing notation would simply be used to achieve that.
Back in Fortran 66 days, you would compile a large 2d array to solve the biggest problem you wanted to solve because you couldn't dynamically allocate. The alternative would be to change the source code and compile, but that was untenable. Then you would set your array elements with ordinary matrix notation but you would need to specify the stride to the solver to work with that: thus yes, it was solving a block, but it is only capable of solving the upper block. So clearly not intended for general submatrix solving: it is meant for using the upper block of a matrix in lieu of dynamically allocating the matrix to the correct size. Alternatively they could have required people to work with 1d arrays and manually handle the matrix indexing, but that would have been unfriendly.
- jabl 5y ago> it was solving a block, but it is only capable of solving the upper block. So clearly not intended for general submatrix solving: it is meant for using the upper block of a matrix in lieu of dynamically allocating the matrix to the correct size. Nope, you can solve any submatrix by passing the address of some other element than the upper left corner as the address of the array.
- BlackFly 5y agoExcept that LOC doesn't take an offset and pointer arithmetic in fortran isn't guranteed to be safe. No practitioner that I ever worked with ever used them besides. I mean you could have simply argued that malloc existed. Anyways, from the original guide of fortran (maybe I am mischaracterizing their rationale, but from the horse's mouth): > However an assumed-size array declaration has been used in the software, in order to overcome some limitations in the Fortran 77 standard. [LAPACK Users' Guide Third Edition] It doesn't matter what the real reason is. The interface to matrix solving algorithms is archaic feeling: DGTTRF( N, DL, D, DU, DU2, IPIV, INFO ) DGTTRS( TRANS, N, NRHS, DL, D, DU, DU2, IPIV, B, LDB, INFO ) Yes, I understand the purpose of every argument, I even understand the weird naming convention. Everything about this feels old and a newer language (such as a newer fortran) could just give you x, _ = A.solve(b) Or something similar and still manage to capture the dynamic sizing, the efficient tridiagaonal storage, sub-block structure if needed, and error information via info; you can trivially make a signature that does it in place if you don't want a silently copying version. There's that LDB, and it certainly has nothing to do with a sub-block of the matrix since it has to do with the storage of the right hand side and enabling a single call to provide multiple solutions while assuming a particular array layout. Anyways, nobody is really engaging with the real point: that it is obvious with use that Fortran is an old language because there are all kinds of these historical artifacts laying around. As I said, everyone knows it, everyone who interacts with it will feel it. But by now fortran isn't alone in this camp. Moreoever, so what? You should be able to admit to yourself that Fortran is an old language (that is simply a fact), that it gives the impression to most people of being an old language (that is again, simply a fact and why shouldn't it?), yet it has many modern features and people can be quite productive in it for some problem domains. Adding modern features isn't going to win people over that are prejudiced against Fortran.
- theodorethomas 5y agoLOC is not Fortran. C_LOC is. "pointer arithmetic" is not Fortran. The Fortran POINTER should not be confused with a C pointer (POINTER has more information). You cannot do arithmetic with a TYPE(C_PTR), either. LDB allows you to solve the 3rd and 5th RHSs only, when B is dimensioned (1000,100): call with NRHS=2, B=B(1,3) and LDB=2000 You can write a type-bound procedure call x = A%solve(b) with today's Fortran. LAPACK calls can be hidden from the casual user. I doubt that a casual user of C++ would enjoy looking at the implementation of its standard library, either. Yes, 29 year-old code written to a 43-year old standard looks old. How is that a Fortran issue? Would it look any newer if it was written in another 40-year old language? The fact that this code is still the mainstay of scientific computation is evidence of the strength of its design and the design of the programming language.
- BlackFly 5y agoFriend, you need to stop missing the forest for the trees.
- jabl 5y agoI meant that instead of passing the array as A or A(1,1) you can pass, say, A(6, 42), if you want to work on the submatrix starting at element (6,42).
- theodorethomas 5y agoWe have moved on, from a compiler that could fit in 4kbytes of main memory to one that is a 2Gigabyte download. It was “array slicing” (array expressions) that was missing in F77. That’s what sequence association was meant to help with. You pass a scalar that is an element of an array and the corresponding dummy is an array. The dummy gets associated with an array that starts at that element and then continues with all the other elements of the actual argument array in array element order. You then use the LDA argument to give a 2d view of the same array in the dummy so you can step to the next column (or every other column if you want, just double the LDA) no matter where you started from. Programs back then had no need for dynamic memory, they were the only thing running on the machine, no kernel, no other processes or users, why not just grab the entire memory? They booted in a blink of an eye, too. Now, not even our supercomputers run like that. You never know who else is using the resources of the machine.
- BlackFly 5y agoYes, I am well aware of what it is doing and we won't really make any progress until you can believe that I know the language. We are disagreeing about the intent and rationale of its addition. Everytime someone tries to interact with DGETRF for the first time, they look up the API (at http://www.math.utah.edu/software/lapack/lapack-d/dgetrf.html http://www.math.utah.edu/software/lapack/lapack-d/dgetrf.htm... maybe) and read, "The leading dimension of the array A. LDA >= max(1,M)." Which tells them nothing, so they keep digging through the documentation. Maybe they check the user's guide, https://www.netlib.org/lapack/lug/node116.html https://www.netlib.org/lapack/lug/node116.html and they read "However an assumed-size array declaration has been used in the software, in order to overcome some limitations in the Fortran 77 standard." but sadly their journey is still not over because none of these source have characterized it as you are (as an array stride due to memory layout). Nor named it appropriately. They still have to keep digging to understand it as a stride. Its name reflects the fact that it was common to oversize the allocation of the memory for the matrix and then LDA > M of course, because in fact LDA = M_MAX. The leading dimension of the matrix in memory (the stride if you were dealing with a 1-d array). If they are lucky (as I was), they will have access to some of the older software that was doing this and they can see what the purpose of the argument and the rationale of the name is. That this also enables one to focus on a sub-block of A and give a stride value. But certainly that isn't the "leading dimension of A". It is the leading dimension of A precisely when A is a bigger matrix and you are focused only on the upper block because of static allocation of A.