title: Fortran memory management and subroutines/functions
body: <p>at the moment I am working on a code for numerical simulations in Fortran 95. My platform is WIndows and I take advantage of the MSVC environment with the Intel Fortran compiler.
  This code, as many in this field, creates a system of equations to be resolved. Numerically, this happen storing a square matrix and a vector of known values. Now, in order to optimize the memory, the matrices are stored in convenient form, like the compressed sparse rows format (CSR) or analogous, so the zero values are not stored.</p>
  
  <p>Given this brief introduction, here there are my doubts.
  Since at compiling time I do not know the dimension of my arrays, I just declare them as:</p>
  
  <pre><code>REAL, DIMENSION(:), ALLOCATABLE :: myArray
  </code></pre>
  
  <p>and once I retrieve the dimension of such a vector, I call</p>
  
  <pre><code>ALLOCATE(myArray(N)) where N is the number of elements that I want to allocate
  </code></pre>
  
  <ul>
  <li>Still, memory is empty, since the values are not stored but a memory check is done in order to avoid stack overflow. Is it right?</li>
  </ul>
  
  <p>Now, filling it with values, the occupied space ramp up. The structure of a Fortran array, both for a 1D vector and multi-dimensional array, is to fill in column order a space equivalent to the number of value. It is to say that if we have a 2D array of dimension 1000x1000, it will be stored in 1M "contiguous boxes" ordered by column numbers (first the first column is stored, then the second one and so on..).</p>
  
  <ul>
  <li><p>If this is true, so the structure of data is the same, is the access time to a particular value the only difference between a multidimensional and a 1D vector?</p></li>
  <li><p>Is then the command <code>RESHAPE</code> changing only the way the program "sees" the arrays?</p></li>
  </ul>
  
  <p>The array I need for my purposes is defined in a module that each subroutine/function share. In particular, a subroutine allocate and fill it. Coming back to the main program, there is no problem with that since I display to the user some statistics about it. Let us say, we allocated 400M REAL*4 values, with about 1.5GB of used memory.</p>
  
  <ul>
  <li><p>However, once that I get into another subroutine, the program stops saying <code>forrtl: severe(170): Program Exception - Stack Overflow</code>. I ran out of memory. But how could it be if the matrix is already allocated and I did not allocate anything more? Notice that: the subroutine uses the same module, so variables are already declared; my RAM has still a free space of about 1.3GB; the stop is at the first line of the subroutine.</p></li>
  <li><p>Is subroutine (and also function) doubling the data? I thought Fortran would pass the address of my variables in that case, avoiding copies and working directly on the values.</p></li>
  </ul>
  
  <p>Finally, as many of you, I enjoyed in C++ the STD library functions, like vector::push_back and so on. In Fortran, there are not such beautiful routines but some very useful functions are still there. Masking an array, using <code>WHERE</code> or <code>COUNT</code> or <code>MERGE</code>can help you to handle some operation effectively.</p>
  
  <ul>
  <li>However, they are veeeeery slow when my matrix is bigger than 1M entries. In that case even a sequential search and substitute is faster than creating a mask or use where. How could it be possible? Aren't they multithreaded?</li>
  </ul>
  
  <p>Thank you in advance for your patience!! All suggestions are very welcome!!</p>