multi dimension dynamic array allocation

I want to use dynamic allocation for multi-block CFD code, where, the index (i,j,k) varies for different blocks. I really dont know, how to allocate arbitrary array index for n number of blocks and pass it to subroutines. I have given a sample code, which gives error message "Error: Expression at (1) must be scalar", on compilation using gfortran.

  common/iteration/nb
  integer, dimension (:),allocatable::nib,njb,nkb
  real, dimension (:,:,:,:),allocatable::x,y,z
  allocate (nib(nb),njb(nb),nkb(nb))
  do l=1,nb
  ni=nib(l)
  nj=njb(l)
  nk=nkb(l)
  allocate (x(l,ni,nj,nk),y(l,ni,nj,nk),z(l,ni,nj,nk))
  enddo
  call gridatt (x,y,z,nib,njb,nkb)
  deallocate(x,y,z,nib,njb,nkb)
  end

  subroutine gridatt (x,y,z,nib,njb,nkb)
  common/iteration/nb
  integer, dimension (nb)::nib,njb,nkb
  real, dimension (nb,nib,njb,nkb)::x,y,z
  do l=1,nb
  read(7,*)nib(l),njb(l),nkb(l)
  read(7,*)(((x(l,i,j,k),i=1,nib(l)),j=1,njb(l)),k=1,nkb(l)),
 $ (((y(l,i,j,k),i=1,nib(l)),j=1,njb(l)),k=1,nkb(l)),
 $ (((z(l,i,j,k),i=1,nib(l)),j=1,njb(l)),k=1,nkb(l))
  enddo
  return
  end

Upvotes: 0

Views: 215

Answers (1)

eriktous
eriktous

Reputation: 6659

The error message gfortran gives, is as good as they get. It points to nib in the line

real, dimension (nb,nib,njb,nkb)::x,y,z

nib is declared as an array. This is not allowed. (What would the size of x, y, and z be in this dimension?)

Apart from this, I don't really understand your description of what it is that you are trying to do, and the sample code you show doesn't make much sense to me.

common/iteration/nb
integer, dimension (:),allocatable::nib,njb,nkb
real, dimension (:,:,:,:),allocatable::x,y,z
allocate (nib(nb),njb(nb),nkb(nb))

When writing new code, using modules to communicate between program units is highly preferred. Old style common blocks are to be avoided.

You are trying to allocate nib, njb, and nkb with size nb. The problem is that nb hasn't been given a value yet (and won't be given one anywhere in the code).

do l=1,nb
ni=nib(l)
nj=njb(l)
nk=nkb(l)
allocate (x(l,ni,nj,nk),y(l,ni,nj,nk),z(l,ni,nj,nk))
enddo

Again the problem with nb not having a value. This loop runs for an unknown number of times. You're also using the arrays nib, njb, and nkb, which don't contain any values yet.

In each iteration of the loop x, y, and z get allocated. This will lead to a runtime error in the second iteration, because you can not allocate an already allocated variable. Even if the allocations would work, this loop would be useless, because the three arrays would be reset in each iteration and would eventually be set to the dimensions of the last allocation.

Now that I'm writing this, I'm starting to think that what it is that you are trying to do, is create so called 'jagged arrays': you want to create a block in x(1,:,:,:) that differs in size in the second, third, and/or fourth dimension from the block in x(2,:,:,:), and so on. This is simply not possible in fortran.

One way to achieve this would be to create a user defined type with an allocatable, three dimensional array component, and create an array of this type. You can then allocate the array component to the desired size for each element of the array of user defined type. This would look something like the following (disclaimer: untested and just one possible way of achieving your goal).

type :: blocktype
    real, dimension(:, :, :), allocatable :: x, y, z
end type blocktype

type(blocktype), dimension(nb) :: myblocks

You can then run a loop to allocate x, y, and z to a different size for each array element. This is assuming nb has been set to the required value, and nib, njb, and nkb contain the desired sizes for the different blocks.

do block = 1, nb
    ni = nib(block)
    nj = njb(block)
    nk = nkb(block)
    allocate(myblocks(block)%x(ni, nj, nk))
    allocate(myblocks(block)%y(ni, nj, nk))
    allocate(myblocks(block)%z(ni, nj, nk))
enddo

If you want to do it like this, you will definitely want to put your procedures in modules, because that way you automatically get explicit interfaces, which are required for passing around such arrays of user defined type.

One afterthought: Don't use implicit typing, not even in sample code. Always use implicit none.

Upvotes: 3

Related Questions